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Elliptical  acoustic  particle  motion  in  underwater  waveguides 

David  R.  Dall’Osto  and  Peter  H.  Dahl 

Mechanical  Engineering  and  Applied  Physics  Laboratory,  University  of  Washington 
1013  NE  40th  St.  Seattle,  WA  98105 

(Dated:  March  22,  2013) 

The  degree  of  circularity  is  a  vector  quantity  relating  to  the  acoustic  particle  motion  induced  in  an 
acoustic  field,  and  is  formulated  by  measurements  of  the  acoustic  particle  velocity  or  the  complex 
vector  intensity.  It  characterizes  the  ellipticity  of  the  particle  motion  path  that,  for  example,  can 
arise  when  there  is  interference  between  a  source  signal  and  its  reflections.  In  an  ocean  waveguide 
the  spatial  dependence  of  this  vector  property  depends  on  the  water  column  sound-speed,  surface 
conditions  and  bottom  properties,  in  addition  to  the  source/receiver  geometry.  In  waveguides  with 
nearly  horizontal  boundaries,  the  vertical  component  of  the  instantaneous  intensity  can  be  used 
to  provide  an  approximation  to  the  degree  of  circularity.  This  approximation,  applied  to  acoustic 
pressure  measurements  from  two  closely  spaced  hydrophones  made  in  2006  off  the  New  Jersey  coast 
in  waters  80  m  deep,  is  used  to  invert  for  environmental  parameters  of  the  sediment  structure. 


PACS  numbers:  43.30.Yj,43.30.Ma,  43.30.Xm,  43.58.Fm 


I.  INTRODUCTION 

In  shallow  water  it  is  common  to  observe  modula¬ 
tions  in  the  amplitude  and  phase  of  the  pressure  field 
due  to  interference  of  an  acoustic  signal  with  reflec¬ 
tions  from  the  sea-surface  and  sea-bed.  Hydrophones, 
which  provide  measurements  of  the  acoustic  pressure 
field,  are  commonly  used  to  interpret  these  modula¬ 
tions  using  pressure-based  processing  schemes  for  source- 
localization,  such  as  beam-forming  with  an  array  of  hy¬ 
drophones,  and  geo-acoustic  inversion.  Vector  sensor 
measurements,  which  include  the  acoustic  pressure  and 
the  three  orthogonal  components  of  the  acoustic  particle 
velocity,  have  been  utilized  in  pressure-based  processing 
schemes  for  source  localization1  as  well  as  geo-acoustic 
characterization.2  For  example,  the  magnitude  and  direc¬ 
tion  of  the  particle  velocity  vector  can  be  readily  trans¬ 
lated  to  the  spatial  dependence  of  the  pressure  phase  for 
use  in  beam- forming  algorithms.3 

Additional  vector  information  can  be  extracted  from 
the  phase  relationship  between  the  individual  compo¬ 
nents  of  particle  velocity;  this  phase  relationship  de¬ 
scribes  a  fundamental  property  of  the  fluid  particle  mo¬ 
tion,  referred  to  as  the  degree  of  circularity.  The  degree  of 
circularity,  which  can  be  formulated  from  measurements 
of  the  acoustic  particle  velocity  or  the  complex  intensity, 
is  a  non-dimensional  number  that  takes  on  characteristic 
values  between  -1  and  1  when  a  signal  interferes  with  its 
reflections.4  The  direction  of  the  particle  velocity  vector 
does  not  necessarily  fall  on  a  straight  line  (equivalent  to 
the  zero  value  for  the  degree  of  circularity)  as  the  parti¬ 
cle  velocity  components  of  interfering  wavefronts  do  not 
necessarily  have  the  same  phase  dependence.  In  ocean 
waveguides,  the  horizontal  surfaces  cause  a  phase  modu¬ 
lation  of  the  vertical  component  of  the  particle  velocity 
with  resulting  elliptical  particle  motion  in  the  vertical 
plane.  This  elliptical  motion  in  the  vertical  plane  can 
be  approximated  by  vertical  line  array  of  closely  spaced 
pressure  sensors. 

We  demonstrate  in  this  paper  how  the  degree  of  circu¬ 


larity  provides  a  useful  metric  to  analyze  the  propagation 
of  a  sound  pulse  in  underwater  waveguides.  The  degree 
of  circularity  has  a  spatial  pattern  that  relates  directly 
to  multi-path  interference.  The  pattern  due  to  the  multi- 
path  interference  depends  on  the  frequency,  acquisition 
geometry,  bottom  properties  and  water  column  proper¬ 
ties.  Spatially  distributed  measurements  of  the  degree  of 
circularity  can  be  used  to  determine  the  depth  and  range 
of  sources,  and  reflections  from  sub-bottom  layers  cause 
a  distinct  pattern  in  range  and  depth  which  can  be  used 
for  geo-acoustic  inversion.  Geo-acoustic  parameters  of  a 
layered  seabed  from  the  Shallow  Water  2006  (SW06)  ex¬ 
periment  near  the  New  Jersey  coastline,  which  has  been 
well  characterized,5  7,9  are  presented  in  light  of  this  vec¬ 
tor  property.  An  inversion  scheme  based  on  the  range 
dependence  of  the  degree  of  circularity  is  presented  along 
with  the  sensitivity  to  specific  geo-acoustic  parameters. 

The  paper  is  organized  as  follows:  Section  II  estab¬ 
lishes  the  mathematical  formulation  for  the  degree  of  cir¬ 
cularity  from  particle  velocity  measurements,  complex 
intensity  measurements,  and  an  approximate  degree  of 
circularity  from  one  directional  component  of  the  com¬ 
plex  intensity.  Section  III  demonstrates  the  effect  of  the 
source-receiver  positions  on  the  degree  of  circularity,  il¬ 
lustrated  by  an  analytical  formulation  of  the  characteris¬ 
tic  Lloyd’s  Mirror  interference  pattern  caused  by  a  source 
below  a  pressure-release  surface.  A  laboratory  experi¬ 
ment  approximating  this  field  demonstrates  the  observ¬ 
ability  of  the  frequency  and  spatial  dependence  of  the 
degree  circularity.  Section  IV  demonstrates  how  the  de¬ 
gree  of  circularity  is  affected  by  reflections  from  layered 
boundaries.  Experimental  data  from  the  SW06  exper¬ 
iment,  from  which  an  approximate  degree  of  circular¬ 
ity  is  formulated  by  vertical  pressure  measurements,  are 
used  to  identify  properties  of  the  sediment  layering  struc¬ 
ture.  Parabolic  equation  (PE)  simulations  are  presented 
to  demonstrate  the  spatial  properties  of  the  degree  of  cir¬ 
cularity.  Section  V  presents  a  parameterized  geo-acoustic 
model  of  the  degree  of  circularity  for  the  SW06  site  to 
demonstrate  the  sensitivity  to  environmental  parameters 
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in  terms  of  the  Cramer  Rao  Lower  Bound  (CRLB).  A 
geo- acoustic  inversion  of  the  SW06  site  based  on  experi¬ 
mental  measurements  of  the  degree  of  circularity  is  also 
presented.  Section  VI  provides  a  summary  of  the  results. 


II.  THE  DEGREE  OF  CIRCULARITY,  © 


A.  ©  computed  via  particle  velocity,  v 


The  magnitude  and  phase  relationship  among  the  in¬ 
dividual  components  of  the  complex  particle  velocity, 
can  be  used  to  construct  a  vector  quantity,  the  degree  of 
circularity, 


0  =  2 


Im(v  x  v*) 


(i) 


where  x  is  the  cross  product,  the  asterisk  is  the  complex 
conjugate,  and  the  factor  of  2  sets  the  bounds  of  any 
component  of  0  between  =bl.  Each  vector  component 
of  0  quantifies  the  degree  to  which  the  particle  motion, 
projected  on  three  orthogonal  planes,  resembles  a  circle. 

When  any  two  components  of  the  particle  velocity  are 
locked  in  phase,  the  particle  motion  in  the  plane  of  the 
two  corresponding  axes  falls  on  a  straight  line.  The  com¬ 
ponent  of  0  describing  this  motion  will  be  perpendicular 
to  this  plane,  and  for  this  straight-line  motion  it  is  equal 
to  zero  as  a  result  of  the  cross  product.  This  is  exempli¬ 
fied  by  a  single  plane  wave  propagating  in  the  x-y  plane 
where  the  z-component  of  0  equals  0. 

If  two  components  of  particle  velocity  are  not  locked  in 
phase,  the  motion  is  elliptical.  The  component  of  0  de¬ 
scribing  this  elliptical  motion  is  between  0  and  ±1,  where 
purely  circular  motion  is  =bl.  This  can  occur  in  acoustic 
fields  characterized  by  multi-path  interference,  such  as 
that  of  a  plane-wave  and  its  reflection  from  a  boundary. 
The  sign  of  0  determines  the  direction  in  which  a  particle 
moves  about  the  path  with  respect  to  the  fixed  coordi¬ 
nate  axis.  Positive  valued  0  represents  counterclockwise 
motion  and  negative  valued  0  represents  clockwise  mo¬ 
tion.  This  distinction  is  further  clarified  in  Sec.  Ill  by  an 
analytic  and  experimental  example  involving  the  Lloyd 
mirror  interference  effect. 


B.  ©  computed  via  complex  intensity,  Ic 

The  degree  of  circularity,  0,  can  be  derived  from  the 
complex  intensity  vector,  Ic.  This  allows  for  an  approx¬ 
imation  based  on  measurements  of  only  one  component 
of  particle  velocity,  which  is  explained  subsequently.  The 
complex  intensity,  which  describes  the  instantaneous  en¬ 
ergy  flux,  is  formulated  by, 

Ic{t)  =  \p(t)v{ty  =I(t)+iQ{i),  (2) 

where  p  is  the  acoustic  pressure,  and  the  real  and  imag- 
inaryparts  of  Ic  represent  represent  the  time  varying 
amplitude  of  the  active  (J)  and  reactive  (Q)  intensity, 
respectively.10  For  a  steady  state  acoustic  field,  Ic  is 
time-independent . 


An  alternative  formulation  of  0  (equivalent  to  Eq.  1) 
based  on  the  complex  intensity  is, 

e-2  JXQ 

\i\2  +  \Q\ 2’ 

Note  that  reactive  intensity  is  very  closely  linked  to  el¬ 
liptical  particle  motion:  when  \Q\  =  0,  the  degree  of 
circularity  is  zero;  when  \Q\  >  0,  the  sign  of  the  reactive 
intensity  relates  directly  to  the  particle  motion  polariza¬ 
tion. 


C.  Approximation  of  ©  component  in  the  vertical  plane 

The  sound  fields  discussed  here  are  assumed  to  ex¬ 
hibit  azimuthal  symmetry  (i.e.  the  particle  motion  lies 
only  in  the  vertical  plane).  Only  one  component  of  the 
vector  0  =  [0r,02,0©]  in  cylindrical  coordinates  is 
non-zero,  the  component  perpendicular  to  the  r-z  plane 
(0©).  Hence-forth,  for  azimuthally  symmetric  fields  we 
will  drop  the  vector  notation  and  simply  denote  the  cir¬ 
cularity  as  0. 

An  approximation  to  0  in  the  vertical  plane  can  be 
derived  from  two  vertically  separated  hydrophones.  This 
approximation  relies  on  an  accurate  estimate  of  the  com¬ 
plex  horizontal  intensity.  The  horizontal  active  inten¬ 
sity  can  be  estimated  by  vector  subtraction  of  the  ver¬ 
tical  component  of  active  intensity  from  the  total  ac¬ 
tive  intensity.  The  total  active  intensity  can  be  approx¬ 
imated  directly  from  the  pressure  measurements  as  the 
pressure  magnitude  squared  divided  by  the  characteristic 
impedance,  poc ,  (or  plane  wave  intensity).  The  resulting 
approximate  horizontal  active  intensity  is 


where  Iz  is  the  vertical  intensity,  and  the  tilde  indicates 
an  approximate  of  the  horizontal  intensity,  Ir. 

As  a  word  of  caution,  in  regions  of  strong  destructive 
interference,  the  magnitude  of  the  total  active  intensity 
may  be  underestimated  by  the  plane  wave  intensity.  This 
is  often  encountered  in  shallow  water,  where  the  multiple 
upward  and  downward  propagating  waves  set  up  stand¬ 
ing  waves  in  the  vertical  direction.  In  regions  of  complete 
destruction  the  magnitude  of  the  vertical  component  of 
particle  velocity  can  exceed  the  magnitude  of  pressure  di¬ 
vided  by  poc.  Measurements  in  these  regions  may  have  a 
negative  value  under  the  square  root  in  Eq.  4,  resulting 
in  an  imaginary  component  of  active  intensity,  an  ap¬ 
parently  nonsensical  result.  Although  the  approximation 
breaks  down,  the  ratio  of  the  vertical  intensity  component 
to  the  plane  wave  intensity  serves  as  an  indicator  of  dis¬ 
locations,  or  nulls,  in  the  field.  In  these  cases,  this  ratio 
serves  as  a  field  indicator  which  can  be  used  to  determine 
the  spatial  pattern  of  nulls,  which  in  shallow  water  are 
highly  dependent  on  the  modes  of  the  waveguide.10 

The  approximation  of  total  intensity  used  in  Eq.  4  is 
appropriate  as  long  as  the  interfering  waves  do  not  cause 
complete  destructive  interference,  i.e.  waves  traveling  in 
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FIG.  1.  (a)  3  dB  contour  intervals  of  the  pressure  magnitude,  relative  to  the  level  at  1  m,  due  to  a  2000-Hz  point  source  1  m 
below  the  air- water  interface.  Exaggerated  particle  paths  at  a  range  of  3.17  m  and  10  depths  are  shown,  for  which  ©  ranges 
from  linear  (0  —  0)  to  circular  (0  =  d=l).  (b)  0  as  a  function  of  depth  at  the  range,  3.17  m.  The  dots  correspond  to  the 
discrete  particle  depths  depicted  in  (a). 


the  same  direction  (e.g.  upward  or  downward).  If  up¬ 
ward  and  downward  traveling  waves  interference,  they 
can  completely  annilate  eachother  and  cause  a  disloca¬ 
tion  where  |p|/|v|  ^  poc.  Furthermore,  the  horizontal 
component  of  reactive  intensity,  Qr,  is  virtually  zero  a 
few  wavelengths  from  the  source,  i.e.  when  kr  »  1, 
where  k  is  the  acoustic  wavenumber  and  r  is  the  slant 
range  to  the  source.  Substituting  Ir  and  Qr  =  0  for  the 
appropriate  quantities  in  Eq.  3,  the  approximate 


0  = 


Qzlr 


Ir 2 


Qz 


(5) 


captures  the  dominate  features  of  the  true  0,  where 
the  tilde  indicates  an  approximate  measure  of  circular¬ 
ity.  Most  importantly,  0  can  be  formed  from  two  closely 
spaced  (<  A/4)  hydrophones,  extending  this  analysis  of 
particle  motion  to  observations  made  by  a  vertical  line 
array  of  hydrophones.  Note  that  the  accuracy  of  the 
approximate  0,  as  computed  by  two  hydrophones,  is  de¬ 
pendent  not  only  on  the  accuracy  of  the  finite- difference 
estimate  of  particle  velocity,  but  also  the  accuracy  of  the 
assumption  that  |p|/|v|  =  poc. 


III.  EFFECTS  OF  SOURCE  POSITION  AND 
FREQUENCY  ON  0 


at  the  source  depth  and  the  other  (a  virtual  source  rep¬ 
resenting  the  reflection)  equidistant  from,  but  on  the  op¬ 
posite  side  of,  the  pressure  release  surface.  The  pressure 
field  is  given  by 


p—ikripiu>t  p—ikr2piwt 

p{x,  y;  t)  =  Po - Po - ,  (6) 

n  r2 

where  ri  =  a/(^  -  re 0)2  +  (y  —  2/o)2, 
and  r2  =  i/(rc-rco)2  +  (y  +  y  o)2, 

and  where  [x,y\  are  the  horizontal  range  and  depth  of 
the  field  point,  [xo,yo\  is  the  position  of  the  source,  uj 
is  the  angular  frequency,  and  po  represents  the  source 
level  (with  units  pressure  times  distance).  The  reflec¬ 
tion  coefficient,  which  modifies  the  second  term  (virtual 
source),  is  —1  due  to  the  pressure  release  boundary  con¬ 
dition.  This  assumption  is  reasonable  at  the  frequencies 
(>2000  Hz)  and  source  depth  (1  m)  used  in  this  demon¬ 
stration.  (A  point  source  with  kzo  <C  1,  where  zq  is  the 
source  depth,  can  radiate  a  significant  portion  of  energy 
across  the  air- water  interface,  a  phenomenon  referred  to 
as  transparency;11  the  values  of  kzo  considered  here  are 
all  greater  than  1.) 

The  acoustic  particle  velocity  field  is  related  directly 
to  the  pressure  field  by  conservation  of  momentum, 


Consider  an  underwater  sound  source  located  close  to 
the  sea  surface.  The  pressure  field  of  a  tone  source  ex¬ 
hibits  the  characteristic  Lloyd’s  mirror  interference  pat¬ 
tern,  which  has  a  straight  forward  analytical  solution. 
The  pressure  field  is  the  sum  of  two  point  sources,  one 


(7) 

at  po 

where  po  is  the  ambient  density  of  the  medium. 
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FIG.  2.  (a)  Experimental  measurements  of  pressure  magnitude  \p\  (solid  black)  and  particle  velocity  magnitude  |u|  (dashed 
black)  versus  depth  at  2000  Hz  are  along  with  the  analytic  solution  of  \p\  (solid  gray)  and  particle  velocity  magnitude  |u| 
(dashed  gray).  Each  profile  is  normalized  by  the  maximum  value,  (b)  Experimental  measurements  of  the  true  circularity  0 
(circle’s)  and  the  approximate  ©  (x’s)  computed  from  the  vertical  component  of  intensity  are  plotted  over  the  analytic  model 
of  ©  (gray).  Results  are  repeated  in  (c)  and  (d)  at  3000  Hz,  and  in  (e)  and  (f)  at  4000  Hz. 


A.  Analytical  depth  dependence  of  © 

The  interference  pattern  of  Lloyd’s  mirror  is  dependent 
on  both  the  source  depth  and  frequency.  Figure  1  shows 
3  dB  contours  of  the  pressure  magnitude  field  (relative  to 
the  maximum)  for  a  2000  Hz  tone  source  at  depth  1  m. 
The  arrows  in  Fig.  1  show  the  path  (highly  exaggerated) 
that  fluid  particles  at  a  range  3.17  m  and  10  representa¬ 
tive  depths  follow  over  3/4  of  a  cycle.  These  paths  are 
exaggerated  by  the  same  factor  to  preserve  their  relative 
magnitudes  (actual  displacements  are  0(  10-11)  m  for  a 
source  level  of  120  dB  re  1/iPa  @  1  m). 

The  depth  dependence  of  0  is  shown  in  the  inset  of  Fig. 
1.  At  the  depths  0.65  m  and  2.25  m,  the  particle  motion 
(and  thus  particle  velocity)  falls  on  a  line  oriented  with 
the  ridge  of  constructive  interference.  A  similar  straight 
line  path  occurs  at  the  trough  of  destructive  interference 
at  1.33  m,  but  here  the  line  is  oriented  in  the  direction 
perpendicular  to  the  trough  of  destructive  interference. 
Particle  motion  that  does  not  fall  on  a  straight  line  is 
observed  on  either  side  of  the  straight  line  regions.  The 
motion  here  is  best  described  as  elliptical. 

As  a  reminder,  the  acoustic  particle  motion  falls  be¬ 
tween  two  extremes  as  parameterized  by  0,  purely  cir¬ 
cular  motion  when  0  =  ±1  and  motion  falling  along  a 
line  when  0  =  0.  The  sign  of  0  signifies  the  direction  of 
travel,  or  polarization,  of  the  particle  motion.  Note  that 
above  or  below  a  position  with  zero  circularity,  the  polar¬ 
ization  abruptly  changes  from  clockwise  to  counterclock¬ 
wise.  For  example,  at  1.13  m  0  is  negative  (clockwise) 
and  at  1.53  m  it  is  positive  (counterclockwise). 

The  sign-change  on  either  side  of  a  point  with  zero 


circularity  can  be  also  be  thought  of  as  indicating  ei¬ 
ther  prograde  or  retrograde  motion.  Prograde  motion 
occurs  when  the  particle  motion  directly  above  and  below 
a  point  of  zero  circularity  converge  and  travel  in  the  same 
direction  as  the  energy  propagation.  Prograde  motion  oc¬ 
curs  along  crests  of  constructive  interference.  Conversely, 
retrograde  motion  occurs  near  destructive  interference, 
evidenced  by  converging  particle  motion  directly  above 
or  below  travels  in  the  opposite  direction  of  energy  prop¬ 
agation.  The  divergence  (or  convergence)  of  the  reactive 
intensity  can  be  used  to  qualify  either  prograde  (diver¬ 
gent  V  •  Q  >  0)  or  retrograde  (convergent  V  •  Q  <  0) 
motion.  Examples  of  divergent  (and  convergent)  reac¬ 
tive  intensity  in  ocean  waveguides,  and  how  they  relate 
to  dislocations  in  the  field,  can  be  found  in  Dall’Osto,  et 
al.10 

Abrupt  changes  in  the  sign  of  0  are  readily  observ¬ 
able  in  experimental  measurements  as  is  discussed  next. 
Furthermore,  the  depth  spacing  of  these  transitions  (and 
their  frequency  dependence)  has  implications  in  source 
depth  and  ranging.  Although  not  developed  here,  these 
abrupt  transitions  relate  directly  to  the  interference  pat¬ 
tern,  which  can  identify  the  spatial  structure  of  the  in¬ 
terference  used  in  source  depth  inversion  algorithms.8 


B.  Experimentally  measured  depth  dependence  of  © 

These  same  particle  motion  features  predicted  by  the 
analytic  solution  (Fig.  1)  have  been  observed  experimen¬ 
tally  at  the  University  of  Washington  Applied  Physics 
Laboratory’s  acoustic  test  facility  (ATF).  An  approxi- 
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mat  ion  to  Lloyd’s  Mirror  was  generated  by  an  acoustic 
projector  (USRD  standard  type  F37  transducer)  posi¬ 
tioned  at  depth  1  m.  Three  pulses  at  center  frequencies 
2,  3,  and  4  kHz,  each  of  4  ms  duration,  were  recorded  by 
an  epoxy  encased  tri-axial  accelerometer  (model  VHS- 
100  from  Ocean  Applied  Acoustic-Tech,  China). 

The  accelerometer  was  suspended  in  a  frame  by  me¬ 
chanical  springs  to  decouple  the  device  from  the  frame 
and  to  maintain  its  orientation.  Two  pressure-sensing 
hydrophones,  spaced  15  cm  apart,  were  attached  to  the 
frame  directly  above  and  below  the  accelerometer.  The 
frame  was  attached  to  a  vertical  profiler  that  moved  the 
sensor  in  10  cm  depth  increments.  To  recreate  Lloyd’s 
Mirror,  the  received  signal  was  gated  to  include  only  the 
first  1  ms  of  interference  between  the  direct  and  surface 
arrival;  after  1  ms,  additional  unavoidable  reflections  in 
the  ATF  (e.g.  pontoons,  lake  bottom)  complicate  the  in¬ 
terference  pattern.  At  each  depth  increment  0  is  calcu¬ 
lated  by  Eq.  1,  and  0  by  Eq.  5.  To  compute  the  complex 
intensity  for  the  circularity  calculations,  the  pressure  at 
the  center  of  the  accelerometer  was  approximated  by  the 
average  of  the  two  hydrophone  signals. 

Figure  2a,  c,  and  e  show  the  experimentally  measured 
magnitude  of  the  pressure  (solid  black)  and  particle  ve¬ 
locity  (dashed  black)  plotted  relative  to  their  maximum 
recorded  level  at  2,  3,  and  4  kHz,  respectively.  These 
data  follow  the  same  depth  dependence  as  the  analytic 
solution  given  in  Eq.  6  for  the  pressure  (solid  gray)  com¬ 
puted  and  the  particle  velocity  (dashed  gray). 

Figure  2b,  d,  and  f  show  the  measured  true  0  (circle’s) 
and  approximate  0  (x’s)  for  2,  3,  and  4  kHz,  respectively. 
These  measurements  agree  with  the  analytic  model  of  0 
(gray),  as  they  exhibiting  either  a  slowly  changing  circu¬ 
larity  though  maxima  in  the  field  or  an  abrupt  change  in 
sign  at  the  destructive  nodes.  The  depths  of  these  fea¬ 
tures  are  related  directly  to  the  frequency  of  the  signal, 
and  the  depth  and  range  of  the  source. 

The  distinct  advantage  of  the  0  description  over  the 
pressure  or  particle  velocity  magnitude,  is  the  ability 
to  distinguish  the  proximity  and  direction  of  a  max¬ 
ima/node  from  one  measurement  point.  That  is,  on  ei¬ 
ther  side  of  a  null  or  maximum,  the  pressure  field  will 
appear  to  have  the  same  magnitude  while  the  value  of 
0  will  be  of  opposite  sign.  A  change  in  the  sign  of  0 
from  one  to  the  other  suggests  passing  though  a  minimum 
(or  maximum)  while  equivalent  pressure  magnitude  mea¬ 
surements  may  indicate  a  constant  pressure  level.  Addi¬ 
tionally,  since  0  is  non-dimensionalized  with  respect  to 
the  source  level,  analysis  based  on  0  is  not  sensitive  to 
source  level  fluctuations.  This  applies  equally  to  0  de¬ 
rived  completely  from  the  particle  velocity  (Eq.  1)  or 
intensity  measurements  (Eq.  3). 

IV.  EFFECTS  OF  A  LAYERED  SEABED  ON  © 

For  a  source  close  to  the  water- sediment  interface,  a 
similar  Lloyd’s  mirror  pattern  (as  discussed  in  the  pre¬ 
vious  section)  can  be  observed.  Unlike  the  sea-surface 
reflection,  where  little  energy  is  transmitted  through  the 
air-sea  interface,  a  portion  of  the  energy  incident  on  the 


TABLE  I.  Models  of  the  sediment  layer  above  the  R-reflector, 
for  nominal  location  39.0245  N,  73.0377  W.  The  layer  is  pa¬ 
rameterized  by  the  thickness,  Hs ,  and  soundspeed,  cs.  Two 
values  indicate  models  with  multiple  overlaying  sediment  lay¬ 
ers,  the  bottom  layer  (R-reflector)  is  modeled  as  a  half-space 
with  soundspeed  1740  m/s.  Broadband  PE  simulations  are 
used  to  compare  the  various  sediment  models  in  terms  of  the 
characteristic  wavelength  of  the  interference  pattern,  A©.  The 
data  show  A©  ~  16  m. 


Soundspeed  (cs) 

Thickness  ( Hs ) 

A©(m) 

Ref. 

1680  m/s 

1630  m/s 

0.2  m 

21.8  m 

12  m 

Choi,  et  al.6 

1581  m/s 

1602  m/s 

14.5  m 

7.4  m 

19  m 

Jiang,  et  al.6 

1605  m/s 

21  m 

16  m 

Jiang,  et  al.7 

1650  m/s 

half-space 

n/a 

Yang,  et  al.9 

1616  m/s 

16  m 

16  m 

this  work 

water-sediment  interface  is  transmitted  into  the  sediment 
and  can  interact,  under  certain  conditions,  with  sub- 
surficial  layers  to  the  extent  they  exist. 

The  particular  environment  used  in  this  discussion  is 
from  the  Shallow  Water  2006  (SW06)  experiment  off  the 
New  Jersey  Coast.  The  sediment  structure  of  this  site 
has  been  the  focus  of  a  number  of  studies,  a  majority  of 
which  agree  on  a  two  layer  model  where  a  sediment  layer 
(thickness  ^20  m)  overlays  an  strong  impedance  contrast 
layer  (referred  to  as  the  R-reflector).  The  presence  of 
this  layer  and  underlying  reflector  is  directly  visible  in 
chirp-sonar  profiling  data.12  Results  from  a  number  of 
geo-acoustic  inversion  studies  from  this  site,  centered  at 
39.0245°  N,  73.0377°  W,  are  summarized  in  Table  I. 


A.  Experimental  description 

The  data  used  in  this  discussion  was  recorded  on 
a  moored  vertical  line  array  (MORAY)  located  at  the 
aforementioned  coordinates.  A  source  (ITC-2010)  posi¬ 
tioned  within  5  m  of  the  bottom  was  towed  from  the 
stern  of  the  R/V  Knorr  starting  at  range  50  m  from 
the  MORAY  and  proceeding  to  range  300  m.  The  ship 
traveled  at  a  speed  of  0.1  m/s  along  four  different  source- 
receiver  bearings  while  the  source  transmitted  a  sequence 
of  three  pulses  every  12  seconds;  a  2-16kHz  FM  sweep, 
a  2  ms  CW,  and  then  a  3  ms  CW,  each  CW  with  cen¬ 
ter  frequency  1  kHz.  The  hydrophone  pair  (separation 
20  cm)  centered  at  the  depth  49  m  is  used  to  compute 
the  approximate  degree  of  circularity,  0.  The  two,  1  kHz 
cw  pulses  are  utilized  for  the  0  calculations,  while  the 
matched-filtered  FM  sweep  data  provides  information  on 
the  arrival  structure.  The  GPS  coordinates  of  the  source 
were  recorded  for  each  transmission. 

Figure  3  shows  the  geometry  of  the  experiment;  the 
ray-paths  of  the  first  three  arrivals  are  shown  for  a  source- 
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FIG.  3.  The  experimental  geometry  and  the  first  three  arrivals.  The  direct  (solid),  bottom  (dashed),  and  subbottom  (dotted) 
path  from  the  towed  source  (depth  75  m)  to  receiving  array  (depth  49  m)  are  shown  for  a  source  range  of  200  m.  The  sound 
speed  profile  in  the  water  column  varies  from  a  maximum  of  1533  m/s  at  the  surface  to  a  minimum  1485  m/s  at  the  thermocline 
at  20  m  depth,  the  sound  speed  at  the  water-sediment  interface  is  1499  m/s. 


receiver  range  of  200  m.  Although  the  amplitude  of 
the  sub-bottom  reflection  (dotted  line)  is  substantially 
weaker  than  the  direct  (solid  line)  and  bottom  (dashed 
line)  arrivals,  the  interference  caused  by  this  arrival  can 
be  observed  in  the  experimental  data.  The  arrival  times 
of  these  first  three  paths,  extracted  from  the  matched- 
filtered  arrivals  of  the  FM-sweep  pulse,  show  that  the  di¬ 
rect,  bottom  and  sub-bottom  path  arrive  simultaneously 
at  a  range  of  200  m.  In  addition  to  the  matched-filtered 
data,  the  measured  0  exhibits  a  range-dependence  that 
we  show  corresponds  to  the  arrival  angle  of  the  sub¬ 
bottom  path. 

B.  PE  simulations  of  ©  for  the  SW06  site 

Broadband  parabolic  equation  (PE)  simulations  of  the 
SW06  site  (computed  via  the  RAM-PE  numerical  code) 
demonstrate  the  effect  of  the  sub-bottom  layer  on  the 
degree  of  circularity.  The  bottom  model  used  in  the  sim¬ 
ulations  presented  in  Fig.  4  consists  of  2  layers:  the  up¬ 
per  layer  properties  are  cs=1614  m/s,  ps  =  2000  kg/m3, 
Hs= 16  m,  the  lower  R-reflector  half-space  properties  are 
Cb= 1740  m/s,  pij  =  2200  kg/m3  (see  Fig.  3  for  variable 
definitions).  Attenuation  in  the  sediment  layer  was  set  to 
0.05  dB/m/kHz,  and  a  value  of  0.3  dB/m/kHz  was  used 
for  the  absorbing  half-space.5 

The  water  column  sound  speed  profile,  cw(z),  used  in 
the  PE  simulations  is  shown  in  Fig.  3.  Time-series  of  the 
pressure  and  particle  velocity  are  computed  via  Fourier 
synthesis  of  the  experimental  pulse  (shown  in  the  inset 
of  4a) ,  which  in  turn  are  used  to  compute  the  time  series 


of  Ic  and  0  at  each  range  and  depth  point.  The  image 
of  0  in  Fig.  4a  is  constructed  by  computing  the  value  of 
0  at  the  time  of  the  direct  arrival. 

At  the  receiver  depth  of  49  m  (white  horizontal  line 
in  Fig.  4a)  the  direct  and  sub-bottom  arrivals  do  not 
interfere  at  ranges  <  150  m.  The  experimental  pulse 
determines  this  range;  for  example,  were  it  to  have  been 
longer  in  time  (or  processing  bandwidth  decreased),  the 
point  at  which  direct  and  sub-bottom  arrivals  begin  to 
interfere  would  be  less  than  150  m.  Since  the  interference 
at  ranges  less  than  150  m  consist  only  of  the  direct  and 
bottom  arrivals,  the  depth  and  range  dependence  of  0 
is  very  similar  to  that  discussed  in  Sec.  II.  If  there  was 
no  layering  structure  in  the  sea-bed,  this  pattern  would 
extend  until  the  range  where  the  surface  path  begins  to 
interfere  with  the  direct  path.  (Additional  interference 
due  to  the  surface  path  can  be  seen  in  Fig.  4  at  depths 
less  than  5  m.) 

Although  the  interference  pattern  for  range  <  150  m 
contains  little  information  on  the  bottom  properties,  the 
depth  and  range  dependence  of  0  relates  directly  to  the 
source-receiver  geometry.  The  range  dependence  of  0 
is  dependent  on  the  angular  difference  of  the  interfering 
arrivals,  the  smaller  the  angular  difference  of  the  two  ar¬ 
rivals,  the  greater  the  scale  of  the  range  dependence.  The 
range  dependence  of  alternating  signed  0  falls  somewhere 
between  [A/2,  oo).  A  model  of  the  range  dependence  is 
given  in  Sec.  V.  The  greater  the  difference  of  the  in¬ 
terfering  arrival  angles,  the  less  sensitive  the  range  de¬ 
pendence  of  0  is  small  changes  in  the  arrival  angles,  as 
could  be  caused  by  source-depth  fluctuation.  Note,  the 
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range  (m)  range  (m) 

FIG.  4.  (a)  Image  of  broad-band  PE  simulated  value  of  0  occuring  during  the  direct  path  arrival  as  a  function  of  range  and 
depth.  See  text  for  bottom  model.  The  range  corresponding  to  the  experimental  data  is  shown  by  white  horizontal  line  at  the 
depth  of  49  m.  (b)  Plot  of  experimentally  measured  ©  over  three  tracks  (thin  black)  and  their  average  (thick  gray),  (c)  Plot 
of  the  corresponding  model  data  (thick  black)  and  the  averaged  data  (thick  gray)  from  (b). 


bottom-direct-path  interference  pattern  is  very  sensitive 
to  changes  in  source  depth  as  the  angular  difference  be¬ 
tween  these  two  arrivals  approaches  zero  with  increasing 
range. 

C.  Experimental  measurements  of  0  at  the  SW06  site 

As  a  result  of  the  ship  motion  over  the  course  of  the 
experiment,  the  interference  pattern  of  the  direct  and 
bottom  arrival  are  not  consistent  over  the  experimental 
runs.  A  shallower  or  deeper  source-depth  will  expand 
or  contract  the  angular  difference  of  the  direct  and  bot¬ 
tom  arrival.  Since  the  angular  difference  of  these  paths 
are  close  to  zero,  the  fluctuations  in  angular  difference 
are  on  the  order  of  the  angular  difference,  which  drasti¬ 
cally  changes  the  range  of  the  interference  nodes  effec¬ 
tively  smearing  out  this  interference  pattern.  The  inter¬ 
ference  pattern  due  to  the  sub-bottom  arrival,  however, 
is  stable  over  the  experimental  runs.  This  insensitivity 
to  source  depth  fluctuation  is  due  to  the  fact  that  the 
propagation  angle  of  the  sub-bottom  path  asymptotes  to 
the  critical  angle  with  increasing  source-receiver  range, 
while  the  direct  and  bottom  propagation  angles  converge 
to  be  nearly  horizontal.  Thus,  the  fluctuations  in  the 
angular  difference  between  the  sub-bottom  path  and  the 
direct/bottom  path  are  only  a  fraction  of  their  total  an¬ 
gular  difference,  which  remains  near  the  critical  angle 
with  increasing  range. 

The  experimentally  measured  0  at  49  m  depth  was 
computed  from  the  pressure  and  vertical  particle  veloc¬ 
ity  (via  the  finite  difference)  as  outlined  in  Section  II- C. 


Figure  4b  shows  0  as  a  function  range,  corresponding 
to  the  value  of  0  at  the  time  of  the  direct  arrival  for  the 
250,  1  kHz  pulses  that  arrive  over  the  range  span  shown  in 
Fig.  4a.  Results  from  the  individual  ship  tracks  along  the 
bearings  120°,  210°,  and  300°  (thin  black  lines)  show  a 
similar  feature  between  175-230  m  that  remains  through 
averaging  (thick  gray  line).  Figure  4c  shows  the  corre¬ 
sponding  PE  model  of  0  (black  line)  plotted  against  the 
averaged  experimental  data  (gray  line). 

As  observed  in  both  the  experimental  and  modeled 
data,  interference  between  the  first  three  arrivals  causes 
distinct  oscillations  in  0  as  a  function  of  range.  At  the 
ranges  where  the  first  three  paths  interfere,  the  0  range 
dependence  can  be  described  by  a  wavelength,  A©,  de¬ 
fined  as  the  characteristic  range  between  peaks  of  similar 
signed  values  of  0.  The  experimental  data  yields  A©  ~  16 
m  and  the  value  of  A©  computed  from  PE  simulations  of 
candidate  bottom-models  are  shown  in  third  column  of 
Table  I.  As  with  any  measurement,  the  measure  of  A©  is 
subjected  to  some  level  of  uncertainty,  which  effects  the 
confidence  in  a  particular  bottom-model. 

V.  SENSITIVITY  OF  A©  TO  BOTTOM  PARAMETERS 

The  measurement  of  A©  from  two  vertically  spaced  hy¬ 
drophones  has  an  inherent  degree  of  uncertainty,  cr\Q. 
This  is  set  by  the  signal  to  noise  ratio  (SNR)  and  an 
analysis  of  PE  simulated  data  with  additive  Gaussian 
white  noise  with  SNR  of  30  dB,  roughly  corresponding 
to  that  in  effect  for  the  field  data,  puts  cr\Q  ~  1.1  m. 
Similarly,  the  standard  deviation  of  A©  measured  over 
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the  three  range-tracks  that  comprise  the  experimental 
data  is  a \&  1  m. 

A  simple  relationship  can  be  used  to  describe  A©  by 
two  propagation  angles,  which  we  exploit  for  a  discussion 
of  the  sensitivity  of  A©  to  the  bottom  parameters  in  terms 
of  the  Cramer  Rao  Lower  Bound  (CRLB)  framework. 
For  a  pulse  with  acoustic  wavelength  (A  m  cw/f ),  this 
relationship  is 


cos  (fid  —  cos  (j)s 

where  <fid  is  the  propagation  angle  of  the  direct  path  and 
( t>s  is  the  propagation  angle  of  sub-bottom  path. 

For  a  given  bottom  parameterization,  ray  theory  can 
be  used  to  compute  the  arrival  angles,  where  the  propa¬ 
gation  angle  of  the  direct  path  is  readily  determined  by 
the  source-receiver  geometry  and  the  propagation  angle 
of  sub-bottom  path  is  constrained  by  Snell’s  law.  This 
constraint  can  be  used  to  construct  the  geometry  of  the 
sub-bottom  path  which  is  sensitive  to  two  parameters, 
the  layer  thickness  ( Hs )  and  the  layer  soundspeed  ( cs ); 
the  density  and  attenuation  of  the  layer  do  not  effect 
the  propagation  angles.  Note  that  a  more  sophisticated 
ray  tracing  method  may  be  necessary  to  compute  these 
angles  if  the  water-column  has  a  significant  sound-speed 
profile.  For  this  analysis,  all  the  paths  considered  re¬ 
main  below  the  thermocline  (see  Fig.  3)  and  simple 
geometrical-methods  are  adequate. 

Figure  5a  shows  the  value  of  A©  for  various  candidate 
two  layer  models  [cs,i7s],  for  which  the  sediment  be¬ 
low  the  layer  is  considered  a  halfspace  with  sound  speed 
1740  m/s.  A  particular  value  of  A©  is  satisfied  by  a  set 
of  [cs,Hs]  parameterizations  and  additional  information 
must  be  included  to  narrow  this  parameter  space. 

One  piece  of  information,  contained  in  the  data,  is  the 
range  dependence  of  the  interference  pattern.  Figure  5b 
shows  contours  of  range,  for  which  the  sub-bottom  path 
and  direct  path  arrival  times  will  be  the  same.  For  ex¬ 
ample,  the  line  identified  with  200  gives  all  combinations 
of  [cs,Hs]  such  that  there  would  be  simultaneous  arrival 
of  sub-bottom  and  direct  paths  at  range  200  m. 

Another  piece  of  information  is  the  two-way  travel  time 
(TWT)  between  the  water-sediment  interface  and  the 
bottom  of  the  layer  from  chirp  sonar  surveys  which  is 
reported  to  be  TWT~  26  ms.6  This  TWT,  which  fol¬ 
lows  the  relationship  TWT  =  2 Hs/cs,  establishes  a  layer 
thickness  of  21-22  m  and  layer  soundspeed  between  1600- 
1700  m/s.  The  dashed  line  in  Fig.  5b  shows  the  parame¬ 
ter  estimates  for  a  TWT^  26  ms  and  the  curved  solid  line 
shows  the  parameter  estimates  for  a  value  of  A©  =  16. 
Finally,  as  noted  in  Sec  IV  A.,  the  FM  sweep  data  has  be 
matched-filtered  to  determine  the  arrival  times  of  the  di¬ 
rect  and  sub-bottom  path,  which  are  equal  at  range  200 
m. 


A.  CRLB  of  parameter  estimates  based  on  0 

The  uncertainty  of  the  bottom  parameters  based  on  a 
value  of  A©  can  be  derived  from  the  Fisher  information 
matrix, 
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FIG.  5.  a)  Contours  of  A©  for  the  source-receiver  geome¬ 
try  shown  in  Fig.  3.  b)  Contours  of  the  range  where  the 
sub-bottom  and  direct  paths  arrive  coincidently.  The  200- 
m  contour  corresponds  to  the  observation  in  the  data.  The 
solid  curved  line  indicates  the  parameter  set  described  by  the 
A©  =  16  m  contour  and  the  dashed  line  indicates  the  parame¬ 
ter  set  from  the  TWT  chirp  sonar  data,  c)  The  error  ellipses 
representing  the  CRLB  for  an  uncertainty  in  the  A©  estimate, 
a  a©  =  1.1  m.  The  grayed- filled  ellipse  corresponds  to  the  un¬ 
certainty  of  the  parameters  for  a  best  fit,  cs  =  1614  m/s  and 
Hs  —  16  m. 
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where  the  elements  are  the  partial  derivatives  of  the  A© 
with  respect  to  the  parameters.  The  CRLB  derived  vari¬ 
ances  and  co- variances  of  the  parameter  set  are  computed 
from  the  inverse  of  the  F  matrix,  F~x .  The  eigenvalues 
of  F~x  are  the  CRLB  estimates  for  the  parameters  under 
study.  If  off-diagonal  elements  of  F~x  are  non-zero  there 
is  evidently  coupling  between  parameter  estimates. 

The  coupling  between  these  two  parameters  is  shown 
by  the  error  ellipses  in  Fig.  5c  for  a  several  parameter 
combinations  [Hs,cs]  given  <ja©  =  1.15  m.  (The  cross 
of  the  major  and  minor  axes  of  the  ellipses  are  centered 
at  [Hs,cs].)  The  orientation  of  the  error  ellipses,  which 
are  determined  by  the  eigenvectors  of  the  F~x  matrix, 
show  that  estimates  of  the  two  parameters  are  coupled, 
and  the  uncertainty  in  Hs  depends  on  the  uncertainty 
in  cs.  For  a  thick  sediment  layer  with  low  sound  speed 
(upper-left  corner)  the  thickness  dominates  the  uncer¬ 
tainty  due  to  the  propagation  angles  being  dominated 
by  the  geometry,  the  refraction  effect  for  a  small  sound- 
speed  contrast  is  less  significant.  The  refraction  effect 
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becomes  more  significant  as  the  propagation  angles  ap¬ 
proach  the  critical  angle  of  the  sediment- water  interface. 
For  a  larger  sound  speed  contrast  and  a  thin  sediment 
layer  (bottom-right  corner),  the  sound-speed  dominates 
the  uncertainty.  Between  these  two  regions,  the  sensitiv¬ 
ity  of  A©  transitions  from  thickness  dominated  to  sound 
speed  dominated,  where  at  some  point  the  parameters 
are  equally  sensitive.  Note  that  this  change  in  relative 
sensitivity  of  the  two  parameters  has  implications  for 
adaptive  estimation  algorithms,13  specifically  the  order 
in  which  parameters  are  estimated. 


B.  Geo-acoustic  inversion  results  based  on  A© 

A  best-fit  parameterization  can  be  constructed  from 
the  intersection  of  the  A©  contour  with  either  the  200  m 
range  contour  or  the  chirp-sonar  constraint  (see  Fig.  5b.) 
Alternatively,  another  best-fit  parameter  estimate  is  one 
based  purely  on  the  TWT  constraint  and  the  200  m  range 
contour.  Regardless  of  the  means  used  to  arrive  at  the 
best  estimate,  the  best  estimate  lies  somewhere  within 
the  triangle  of  intersecting  lines  (see  shaded  region  in 
Fig.  5b). 

We  conclude  that  the  region  defined  by  the  three  pa¬ 
rameter  sets  is  in  agreement  with  the  various  parameter 
estimates  in  Table  I.  The  best  parameter  estimate  in¬ 
ferred  from  the  towed  source  transmissions  is  Hs  =  16 ±4 
m  and  cs  =  1616T10  m/s.  The  possible  values  of  [cs,  Hs\ 
due  to  the  uncertainty  are  contained  in  the  filled  ellipse 
in  Fig  5c.  This  error  ellipse  encompasses  the  majority  of 
the  parameter  space  given  by  shaded  region  in  Fig.  5b 
and  we  deem  this  a  reasonable  best-fit  estimate. 


VI.  SUMMARY  AND  DISCUSSION 

A  vector  description  of  the  acoustic  field  contains  use¬ 
ful  information  relating  to  the  propagation  of  sound. 
As  specifically  discussed  in  this  work,  the  frequency 
and  depth/range  dependence  of  interference  patterns  can 
identified  by  a  vector  property  which  characterizes  parti¬ 
cle  motion,  the  degree  of  circularity  (0).  Measurements 
of  0  are  useful  in  that  the  range,  depth  and  frequency 
dependence  of  0: 

•  relates  directly  to  source  and  receiver  position,  hav¬ 
ing  utility  in  source  localization 

•  relates  directly  to  reflections  from  sediment  layers, 
having  utility  in  geo-acoustic  inversion 

•  can  be  approximated  in  the  vertical  plane  by  a  ver¬ 
tical  line  array  of  hydrophones,  and  thus  can  be  ap¬ 
plied  to  conventional  pressure-based  measurement 
systems 

Furthermore,  in  addition  to  pressure-based  algorithms 
used  in  localization  and  geo-acoustic  inversion,  the  out¬ 
put  of  vector  sensors  can  be  presented  in  light  of  vector 
properties  of  the  sound  field. 


The  results  of  the  geo-acoustic  inversion  based  on  a  the 
range  dependence  of  0  presented  here  are  in  agreement 
with  previously  reported  pressure-based  inversion  results 
of  the  SW06  site.  Although  this  parameterization  is  con¬ 
sistent  with  the  layered  bottom  models  presented  in  the 
literature,  one  study  does  suggests  a  half-space  model.9 
It  is  worth  noting  that  the  circularity  data  collected  along 
the  same  30°  bearing  could  not  be  distinguished  from  the 
noise.  As  the  inversion  technique  presented  here  depends 
the  spatial  dependence  of  particle  motion  circularity,  the 
lack  of  a  signal  along  this  track  supports  this  conclusion 
of  a  half  space  model  for  this  bearing.  Nevertheless,  the 
data  collected  on  the  other  three  bearings  (radials)  sup¬ 
port  the  a  sub-bottom  reflector  model.  The  observations 
of  circularity  accommodate  both  of  these  models,  and 
suggests  that  the  detection  of  sub-bottom  layers  might  be 
due  to  the  sensitivity  of  the  measurement  system,  such 
as  the  signal  to  noise  ratio. 

The  SNR  of  the  circularity  is  dependent  on  the  back¬ 
ground  noise  of  both  the  pressure  field  and  the  particle 
velocity  field.  The  uncertainty  in  the  0  measurements 
could  be  much  improved  upon  by  changing  the  sensing 
technique,  as  the  approximate  technique  used  here  is  sub¬ 
ject  to  the  SNR  of  the  pressure  field.  The  background 
noise  measured  by  an  accelerometer  is  inherently  differ¬ 
ent  than  the  noise  on  a  hydrophone.  Depending  on  the 
relative  SNR  of  the  pressure  or  particle  velocity  field,  the 
calculation  of  0  via  a  particle  motion  sensor  may  be  re¬ 
sult  in  a  higher  SNR,  and  thus  lower  uncertainty  in  the 
characteristic  scales  of  0.  Under  the  same  environmental 
conditions,  this  improvement  in  SNR  will  reduce  the  un¬ 
certainty  in  parameter  estimates  over  those  derived  from 
pressure-based  formulations  of  circularity. 
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Acoustic  intensity  is  a  vector  quantity  described  by  collocated  measurements  of  acoustic  pressure 
and  particle  velocity.  In  an  ocean  waveguide,  the  interaction  among  multipath  arrivals  of  propagating 
wavefronts  manifests  unique  behavior  in  the  acoustic  intensity.  The  instantaneous  intensity,  or  energy 
flux,  contains  two  components:  a  propagating  and  non-propagating  energy  flux.  The  instantaneous 
intensity  is  described  by  the  time-dependent  complex  intensity,  where  the  propagating  and 
non-propagating  energy  fluxes  are  modulated  by  the  active  and  reactive  intensity  envelopes, 
respectively.  Properties  of  complex  intensity  are  observed  in  data  collected  on  a  vertical  line  array 
during  the  transverse  acoustic  variability  experiment  (TAVEX)  that  took  place  in  August  of  2008, 

17  km  northeast  of  the  Ieodo  ocean  research  station  in  the  East  China  Sea,  63  m  depth.  Parabolic 
equation  (PE)  simulations  of  the  TAVEX  waveguide  supplement  the  experimental  data  set  and  pro¬ 
vide  a  detailed  analysis  of  the  spatial  structure  of  the  complex  intensity.  A  normalized  intensity  quan¬ 
tity,  the  pressure-intensity  index,  is  used  to  describe  features  of  the  complex  intensity  which  have  a 
functional  relationship  between  range  and  frequency,  related  to  the  waveguide  invariant.  The 
waveguide  invariant  is  used  to  describe  the  spatial  structure  of  intensity  in  the  TAVEX  waveguide 
using  data  taken  at  discrete  ranges.  ©  2012  Acoustical  Society  of  America. 

[DOI:  10.1121/1.3682063] 

PACS  number(s):  43.30.Yj,  43.30.Bp,  43.58.Fm  [AJZ]  Pages:  2023-2035 


I.  INTRODUCTION 

The  measurement  of  an  acoustic  intensity  vector  field  is 
a  standard  technique  used  for  localizing  sources  and  deter¬ 
mining  sound  power.1  For  example,  active  intensity  is  the 
net  acoustic  power  per  unit  area;  measurements  of  the  com¬ 
ponent  normal  to  a  control  surface  surrounding  a  source 
upon  integration  yields  the  radiated  power.  The  scaled  range 
kr  1,  where  k  is  the  acoustic  wavenumber  and  r  is  range 
from  source,  provides  one  kind  of  far  field  definition  which 
upon  being  satisfied  the  radiation  impedance  becomes  purely 
real  and  intensity  is  considered  active.  This  is  most  readily 
seen  in  the  case  of  a  point  source  where  the  radiation  imped¬ 
ance  equals  p0c(l  +  i/kr),  where  p0  is  the  density  and  c  is 
the  sound  speed  of  the  medium. 

In  a  waveguide,  just  as  in  room  acoustics2  and  other 
non-free  field  environments  involving  multipath  propaga¬ 
tion,  a  scaled  range  kr  1  does  not  necessarily  define  a 
range  beyond  which  the  radiation  impedance  is  purely  real 
and  this  has  implications  on  the  intensity  vector  field.  In  this 
paper,  properties  of  the  acoustic  intensity  vector  field  in  a 
shallow  water  waveguide  are  studied,  using  a  combination 
of  field  data  obtained  from  an  experiment  staged  off  the 
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coast  of  Korea  within  the  East  China  Sea,  and  synthetic  data 
generated  by  parabolic  equation  (PE)  simulations. 

The  paper  is  organized  as  follows.  Section  II  establishes 
the  definitions  (and  notation)  for  instantaneous  and  complex 
intensity  to  be  used  in  this  paper,  including  the  real  (active) 
and  imaginary  (reactive)  parts  of  complex  intensity.  Section  III 
describes  the  field  experiment  and  presents  examples  of  both 
instantaneous  and  complex  intensities  estimated  from  a  vertical 
line  array  used  in  the  field  experiment.  This  section  ends  with 
a  discussion  on  processing  bandwidth  and  its  relation  to 
estimating  complex  intensity. 

In  Sec.  IV  the  focus  turns  to  narrow  band  intensity,  or  in¬ 
tensity  estimated  in  the  frequency  domain.  Following  the  work 
of  Eliseevnin  et  al .,3  a  Pekeris  waveguide  model  is  used  to 
obtain  some  initial  insights  into  the  vector  intensity  properties 
of  ocean  waveguides,  where  this  paper  extends  their  analytic 
analysis  to  normalized  forms  of  the  complex  acoustic  intensity 
vector.  In  this  section  experimental  data  is  again  examined  but 
with  sufficiently  narrow  bandwidth  such  that  intensity  proper¬ 
ties  lose  their  time  dependence.  Using  a  measured  sound  speed 
profile  (SSP)  and  an  empirically  derived  geoacoustic  model 
for  the  seabed,  equivalent  synthetic  results  are  generated  with 
PE  simulation;  properties  of  the  active  and  reactive  intensities 
observed  in  both  the  field  and  synthetic  data  are  discussed  and 
interpreted.  Part  of  this  interpretation  includes  utilization  of 
the  waveguide  invariant,  as  demonstrated  for  reactive  intensity 
in  D’ Spain  et  al.,4  for  particle  velocity  and  active  intensity  in 
Ren  et  al.,5  and  demonstrated  in  this  paper  for  both  active  and 
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(5) 


reactive  normalized  intensity  quantities.  A  summary  is  pro¬ 
vided  in  Sec.  V. 

II.  ACOUSTIC  INTENSITY 

Instantaneous  acoustic  intensity,  S(t),  also  referred  to  as 
the  Umov-Poynting  vector,6  describes  the  instantaneous 
flow  of  energy  in  an  acoustic  field.  To  estimate  5(0  requires 
measurement  of  both  the  pressure  and  particle  motion  at  a 
field  point.  This  quantity  can  easily  be  formulated  from  the 
output  of  a  vector  sensor,  a  combined  receiver  measuring  the 
acoustic  particle  velocity  and  pressure.  The  acoustic  particle 
velocity  can  also  be  formulated  by  a  spatial  distribution  of 
closely  spaced  pressure  sensors,  with  respect  to  the  signal 
wavelength.7  The  equivalence  of  these  measurement  systems 
relies  on  the  fact  that  in  Euler’s  equation 

Vp(0  =  —p0a(t),  (1) 


V(t)  =  VM(t)  +  !%(/), 

where  the  carrot  represents  the  Hilbert  transformed  time- 
series.  These  can  in  turn  be  used  to  construct  the  complex 
intensity, 

!c{t)  =  ^P(t)v*{t),  (6) 

where  the  asterisk  in  this  equation  denotes  the  complex  con¬ 
jugate.  The  complex  intensity,  as  formulated  by  Heyser,8  and 
as  used  in  Eq.  (6),  is  a  time  dependent  quantity.  This  differs 
from  a  frequency  transformed  formulation  of  complex  inten¬ 
sity,  which  is  a  time-independent  definition  often  described 
in  acoustic  intensity  literature.  The  real  part  of  the  complex 
intensity  represents  the  envelope  of  the  active  intensity 

7(0=Re{/c(0},  (7) 


where  a  is  the  particle  acceleration  and  p  is  the  acoustic  pres¬ 
sure,  the  pressure  gradient  can  be  approximated  along  the 
axis  of  two  pressure  sensors  using  the  finite  difference  of  the 
pressure  measurements.  The  particle  velocity,  v,  is  obtained 
by  time  integration 


v{t)  =  - 


rm*'. 

0  Po 


(2) 


The  lower  limit  in  this  case  represents  the  onset  of  the  pres¬ 
sure  signal. 

The  instantaneous  acoustic  intensity  is  defined  by  the 
product  of  real-valued  acoustic  pressure  and  particle  veloc¬ 
ity.  For  complex  representations,  5(0  is  defined  at  a  field 
point  using  the  real  parts  of  these  quantities  as  in 


5(0  =Re{/?(0}Re{v(f)}. 


(3) 


Note  that  complex  representations  p(t)  and  v(t)  can  be  gener¬ 
ated  from  their  measured  real- valued  counterparts,  pM( 0  and 
vM(t ),  by  using  the  Hilbert  transform  of  the  real-valued  sig¬ 
nal  as  discussed  further  below. 

The  instantaneous  energy  flux  through  a  fluid  parcel 
consists  of  the  sum  of  two  components,  active  and  reactive 
intensity.  Active  intensity  is  a  measure  of  the  energy  flux 
propagating  through  a  volume  and  corresponds  to  the  degree 
to  which  the  pressure  and  particle  velocity  signals  are  in- 
phase.  The  time-average  of  the  active  intensity  represents 
the  mean  energy  flow.  Reactive  intensity  is  the  energy  flux 
which  in  a  time- average  sense  is  zero.  The  existence  of  reac¬ 
tive  intensity  depends  on  the  degree  to  which  the  pressure 
and  particle  velocity  signals  are  out  of  phase. 

The  phase  properties  of  an  acoustic  signal  are  effec¬ 
tively  described  using  the  Hilbert-transform  pair  to  form  an 
analytic  signal,  consisting  of  the  original  time-series  [pM( t) 
and  vM(t)]  and  an  imaginary  component  that  is  90°  out  of 
phase,  or  phase-quadrature.  The  analytic  pressure  signal  is 

p{t)  =  PM{t)  +  ipM(t)  (4) 

and  particle  velocity  signal  is 


and  the  imaginary  part  represents  the  envelope  of  reactive 
intensity 

Q(t)  =  lm{/c(0}.  (8) 

These  intensity  envelopes  represent  the  low-pass  filtered 
active  and  reactive  intensity  components.9  Provided  the 
bandwidth  B  is  sufficiently  narrow  with  respect  to  the  center 
frequency,  /0,  such  that  B/fo  <C  1,  then  5(0  can  be  repre¬ 
sented  by  the  sum 

5(0  =  /(0  cos2[2nf0t]  +  Q(0  sin[4nfot\/2.  (9) 

Properties  of  5(0  described  by  Eq.  (9)  are,  for  example, 
readily  observed  in  the  instantaneous  axial  intensity  originat¬ 
ing  from  a  speaker10  in  free  field  conditions;  for  kr  1  nar¬ 
row  band  measurements  (here  third-octave  band,  centered  at 
/o  =  250Hz)  go  as  cos2[27i/00  modulated  by  7(0,  which 
traces  time  variation  in  the  source  power.  The  acoustic  inten¬ 
sity  field  in  this  region  is  active  and  the  time-average  (5(0) 
is  non-zero.  For  kr  <C  1  the  instantaneous  intensity  goes  as 
sin[47i/o0/2  modulated  by  Q(t),  and  here  the  field  is  reactive 
and  (5(0)  is  nearly  zero.  For  intermediate  kr,  the  intensity 
behavior  is  a  combination  of  these  two  fields.  While  this 
behavior  as  a  function  of  kr  exists  in  free  field  conditions, 
multipath  interference  in  a  waveguide  represents  a  source  of 
reactive  intensity  that  exists  even  for  large  kr. 

III.  EXPERIMENTAL  DATA  SET 

The  experiment  was  conducted  17  km  northeast  of  the 
Korean  ocean  research  station  (Ieodo)  in  waters  63  m  deep, 
and  formed  a  component  of  a  larger  program  known  as  the 
transverse  acoustic  variability  experiment  (TAVEX).  Acous¬ 
tic  pressure  data  were  collected  on  a  16-element  vertical  line 
array  (VLA)  with  uniform  hydrophone  spacing  (d)  of  3  m. 
Each  hydrophone  was  sampled  at  a  rate  of  20.47  kHz.  The 
source  signals  were  implosions  of  light  bulbs  (air  volume 
700  cm3)  at  a  nominal  depth  of  40  m  off  the  stern  of  the  R/V 
Sunjin.  Due  to  wind  coupling  of  the  research  vessel  and  the 
resulting  relative  motion  of  the  vessel  with  respect  to  the 
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water  column,  the  bulb  deployment  system  is  suspended  at 
an  angle  from  vertical  resulting  in  bulbs  implosion  depths 
estimated  to  be  37  m.  The  key  feature  of  the  implosion  signal 
is  a  long-time  scale  rarefaction  signal,  followed  by  a 
impulsive-like  spike  in  pressure.  This  sequence  (bubble 
pulse)  is  repeated  every  4.7  ms  with  an  exponential  decay, 
similar  to  that  described  by  Heard. 1 1 

Figure  1(a)  shows  the  locations  of  the  bulb  sources  rela¬ 
tive  to  the  moored  VLA.  The  arrows  represent  the  depth- 
average  current  vector  during  each  bulb  implosion,  inferred 
from  a  two  week  acoustic  doppler  current  profiler  (ADCP) 
record  collected  within  1  km  of  the  VLA  site.  Although  the 
current  measurements  began  the  day  after  the  acoustic  meas¬ 
urements,  the  daily  tidal  signal  in  the  ADCP  data  shifts  very 
little  over  the  two  week  measurement  period,  and  thus  the 
current  properties  during  the  two  hour  period  of  the  acoustic 
measurements  on  August  21,  2008  are  inferred. 

During  the  time  period  of  the  experiment,  the  magnitude 
of  the  current  is  ^0.6  m/s  and  slowly  shifts  direction  from 
NW  to  N.  The  VLA  is  subjected  to  a  depth  dependent  tilt  in 
the  direction  of  the  current.  Estimates  of  the  relative  receiver 
locations  from  ray  arrival-times  of  the  four  closest  sources 
suggest  that  the  array  shape  remains  relatively  constant  with 
respect  to  the  current  direction  [Fig.  1(b)].  Environmental 
data  from  the  nearby  Ieodo  tower  indicates  a  7  m/s  southerly 
wind,  and  seas  characterized  by  an  rms  wave-height  (Hrms) 
of  0.475  m.  SSPs  measured  from  three  expendable  bathyther¬ 
mograph  (XBT)  casts  during  the  experiment  show  a  down¬ 
ward  refracting  environment. 


A.  Instantaneous  and  complex  vertical  intensity 

A  reliable  estimate  of  the  instantaneous  intensity  based 
on  the  finite  difference  approximation  requires  low-pass  fil¬ 
tering  of  the  pressure  time-series  before  computing  the  parti¬ 
cle  velocity  to  satisfy  d/A  requirements.  To  demonstrate  the 
behavior  of  the  instantaneous  intensity  for  the  first  three 
arrivals,  the  filter  bandwidth  must  be  wide  enough  to  resolve 
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FIG.  1.  (a)  Source  locations  (circles)  relative  to  the  VLA  (square)  plotted 
with  their  estimated  ocean  current  vector.  The  labels  A1-A12  denote  the 
time  sequence  of  the  bulb  implosions,  occurring  over  a  two  hour  period,  (b) 
Estimated  array  shape  (relative  hydrophone  offset)  corresponding  to  a  tilt 
directed  along  the  current  vector  for  each  bulb  source 


the  time- separated  arrivals.  Typically,  finite  difference  con¬ 
straints  on  frequency  for  a  given  sensor  spacing  are  set  to 
meet  an  acceptable  percent  error  in  an  unknown  acoustic 
field. 

The  error  for  the  finite  difference  approximation  is  a 
function  of  frequency  and  arrival  angle.  In  situations  where 
the  field  can  be  described  as  locally  planar,  the  finite  differ¬ 
ence  error  in  particle  velocity  is 


2  sin 

kd  sin  9 

2 

kdsinO 

(10) 


and  for  pressure  is 


e(p)  =  1  —  cos 


kdsinO 

2 


(11) 


The  error  in  the  intensity  component  along  the  axis  separat¬ 
ing  the  two  sensors  is  the  sum  of  Eqs.  (10)  and  (11),  and 
only  effects  the  magnitude  of  intensity.9  There  is  no  inherent 
phase  error  in  estimating  pressure  and  vertical  velocity  and 
the  existence  of  reactive  intensity  should  not  be  attributed  to 
errors  in  the  finite  difference  technique.  Measurements  may 
suffer  from  other  sources  of  error  due  to  imperfections  in  the 
measurement  system,  such  as  error  in  the  active  intensity 
due  to  phase  mismatch  between  two  sensors  or  error  in  the 
reactive  intensity  due  to  amplitude  mismatch  between  two 
sensors.12  Fortunately,  the  errors  discussed  here  are  not  an 
issue  when  considering  simulated  data,  which  in  turn  vali¬ 
date  observations  in  the  experimental  data. 

By  limiting  the  analysis  to  propagation  angles  near 
broadside  (where  broadside  represents  9  =  0°),  the  finite  dif¬ 
ference  error  can  be  kept  below  the  aliasing  frequency,  pre¬ 
serving  the  correct  sign  of  pressure  and  particle  velocity. 
Here,  the  aliasing  frequency  occurs  when  the  sign  of  the 
pressure  estimate  at  the  midpoint  used  in  the  finite  difference 
approximation  is  no  longer  preserved,  most  easily  inter¬ 
preted  for  an  endhre  arrival  (0  =  90°)  when  d  =  A/2.  Consid¬ 
ering  wavelengths  that  are  small  in  comparison  to  the 
waveguide  depth,  ray  theory  can  be  employed  to  determine 
characteristics  of  the  first  three  arrivals.  Arrival  angles  of  the 
first  three  eigenrays  for  the  closest  source  location  (Al), 
based  on  the  experimental  SSP  in  Fig.  2(a)  and  shown  in 
Fig.  2(b),  correspond  to  an  average  finite  difference  error  in 
intensity  at  360  Hz  of  2%  for  the  direct  path  (black),  22%  for 
the  surface  path  (light  gray),  and  41%  for  the  bottom  path 
(dashed),  all  below  the  aliasing  frequency.  Note  that  this 
source  range  is  particularly  interesting  in  that  it  provides  a 
time-series  that  exhibits  temporal  separation  between  the 
peak  arrivals  of  the  first  three  acoustic  paths  as  well  as  a 
time-series  that  exhibits  multi-path  interference. 

A  simulated  pressure  times-series  of  the  VLA  data  for 
the  Al  source  is  constructed  to  evaluate  the  effect  of  the 
downward  refracting  SSP  and  multi-path  interference  on  the 
received  signal.  The  simulated  time- series,  B  =  300  Hz  and 
/o  =  210Hz,  is  generated  from  a  Fourier  synthesis  of  para¬ 
bolic  equation  (PE)  solutions  (these  methods  are  further 
described  in  Sec.  IV).  Figure  2(c)  shows  contours  of  \p\2  as  a 
function  of  time  and  depth  receiver  depth.  The  source/ 
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FIG.  2.  (a)  The  average  SSP  (black)  used  in  numerical  simulations  is  derived  from  three  experimentally  measured  SSPs  (gray),  (b)  Eigenrays  of  the  direct 
path  (solid  black),  surface  path  (gray),  and  bottom  path  (dashed)  for  the  midpoints  of  the  15  hydrophone  pairs  (circles),  (c)  Contours  of  \p\ 2  of  a  simulated  bulb 
pulse,  B  =  300  Hz  and /0  =  210  Hz,  formed  by  a  Fourier  synthesis  of  parabolic  equation  (PE)  numerical  solutions.  Arrival  times  of  the  first  three  eigenrays  (dot¬ 
ted),  computed  from  ray  theory,  are  plotted  over  the  contours. 


receiver  ranges  of  the  VLA  sensors  are  dependent  on  the 
array  tilt;  the  approximate  tilt  during  the  A1  signal  is  incor¬ 
porated  into  this  simulation.  The  arrival  times  of  the  first 
three  eigenrays  are  plotted  (dotted  line)  over  the  contours  to 
highlight  the  arrival  times  of  these  paths.  There  are  two 
kinds  of  gradients  in  \p\ 2,  beyond  that  accounted  for  by 
spherical  spreading,  observed  along  the  vertical  axis  of  the 
simulated  times  series.  One  is  due  to  constructive  and  de¬ 
structive  interference  due  to  multi-path  interference;  the 
other  is  due  to  refraction,  owing  to  the  SSP,  which  generates 
a  focusing  and  defocussing  of  the  wavefronts.  As  reactive  in¬ 
tensity  is  directed  along  decreasing  \p\ 2  contours,4  we  antici¬ 
pate  a  degree  of  reactive  intensity  will  be  observable  in  these 
two  regions. 

With  results  of  Fig.  2  in  mind,  we  now  examine  the  field 
data.  For  each  pair  of  neighboring  receivers,  the  instantane¬ 
ous  and  complex  intensity  component  directed  along  the 
array  are  formed  by  Eqs.  (3)  and  (6).  In  the  coordinate  sys¬ 
tem  used  here,  placing  the  origin  at  the  top  of  the  array 
(14.5  m),  negative  intensity  indicates  propagating  energy 
flux  toward  the  surface.  Figures  3(a)  and  3(b)  show  a  nor¬ 
malized  component  of  S{t)  (dashed  lines)  computed  from  the 
hydrophone  pair  centered  at  32.5  m  depth.  The  hydrophone 
signals  were  first  processed  by  a  bandpass  filter  ( B  =  300  Hz 
and  /o  =  210  Hz),  which  were  in  turn  used  to  compute  the 
intensity  component  along  the  array.  These  intensity  signals 
are  normalized  by  the  maximum  squared  pressure  divided 
by  p0c. 

At  this  depth,  the  component  of  S(t )  along  the  array 
starts  off  as  being  predominately  negative,  as  the  energy  flux 
associated  with  the  direct  path  passes  the  sensors.  As  the  sur- 
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face  reflected  and  bottom  reflected  paths  arrive  at  the  re¬ 
ceiver  pair  [nearly  simultaneously  according  to  the  ray 
arrival  times  in  Fig.  2(c)],  the  instantaneous  intensity 
becomes  oscillatory,  i.e.,  with  reduced  time  average,  which 
is  a  characteristic  of  interference  between  the  two  paths.  If 
the  surface  and  bottom  arrivals  were  sufficiently  time  sepa¬ 
rated,  the  instantaneous  intensity  would  appear  to  be  on  one 
side  of  the  axis  (positive  for  the  surface  arrival  and  negative 
for  the  bottom  arrival). 

The  active  intensity,  normalized  in  the  same  manner  as 
the  instantaneous  intensity,  is  plotted  in  Fig.  3(a)  over  the  in¬ 
stantaneous  intensity  by  the  thick  black  line.  Similarly  the 
reactive  intensity  envelope  is  plotted  in  Fig.  3(b)  as  the  thick 
gray  line.  The  negative  active  intensity  indicates  energy  flux 
traveling  up  the  array  during  the  direct  arrival,  which  then 
switches  to  positive  indicating  energy  flux  traveling  down 
the  array  (away  from  the  surface)  during  the  subsequent  ar¬ 
rival.  Even  though  there  is  energy  traveling  upwards  from 
the  bottom  arrival,  the  signal  is  dominated  by  the  surface 
path.  The  interference  between  these  opposing  energy  fluxes 
produces  a  strong  component  of  reactive  intensity. 

There  also  is  a  small  component  of  reactive  intensity 
associated  with  the  direct  arrival  observed  in  the  data  at  5  ms 
in  Fig.  3(b).  This  component  of  reactive  intensity  may  be 
due  to  refraction  and  subsequent  focusing  associated  with 
the  SSP,  as  seen  in  contours  of  \p\ 2  in  the  PE  simulation  in 
Fig.  2(c).  The  SSP-dependent  focusing,  which  has  the  great¬ 
est  effect  on  the  direct  path,  does  not  produce  the  level  of 
reactive  intensity  as  observed  during  interference  between 
the  surface  and  bottom  paths.  The  sign  of  the  reactive  inten¬ 
sity  envelope  does  not  necessarily  describe  the  direction  of 
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FIG.  3.  (a)  The  component  of  instantaneous  intensity  (dashed)  and  the  active  intensity  envelope  (solid  black)  oriented  along  the  array  at  depth  32.5  m  for  the 
bulb-pulse  at  190  m  range,  (b)  The  same  component  of  the  instantaneous  intensity  (dashed)  plotted  over  the  reactive  intensity  envelope  (solid  gray).  The  inten¬ 
sities  are  normalized  by  the  maximum  value  of  \p\2/p0c.  (c)  Image  of  the  active  intensity  envelope  measured  along  the  length  of  the  array,  (d)  Image  of  the  re¬ 
active  intensity  envelope  measured  along  the  length  of  the  array.  Whited  out  regions  correspond  to  a  6  m  gap  in  data  due  to  an  inoperable  hydrophone,  (e) 
Image  of  the  active  intensity  envelope  along  the  modeled  VLA  array  shape  for  a  PE  simulated  bulb  pulse,  (f)  Image  of  the  corresponding  reactive  intensity  en¬ 
velope.  Active  and  reactive  intensity  images  are  normalized  by  the  depth  average  of  the  maximum  value  of  \p\2/po  c. 


the  instantaneous  energy  flux,  but  it  does  provide  an  indica¬ 
tion  of  the  proximity  and  direction  of  gradients  in  \p\2. 

These  intensity  features,  created  from  both  focusing  and 
multi-path  interference,  have  a  spatial  structure  observable 
along  the  length  of  the  VLA.  Figures  3(c)  and  3(d)  show  the 
component  of  the  active  and  reactive  intensity  envelopes 
measured  along  the  array,  respectively,  for  the  A1  bulb 
range.  In  this  case,  the  intensity  envelopes  are  normalized  by 
a  constant  equal  to  the  depth  average  of  the  maximum  mag¬ 
nitude  squared  pressure  divided  by  p0  c •  As  mentioned  ear¬ 
lier,  the  active  intensity  envelope  denotes  the  direction  of 
propagating  energy  flux,  negative  indicating  energy  flux  up 
the  array  toward  the  surface.  Figures  3(e)  and  3(f)  show  the 
active  and  reactive  intensity  envelopes  of  a  simulated  bulb- 
pulse  measured  along  the  tilted  array.  The  simulated  bulb- 
pulse  is  from  a  Fourier  synthesis  of  PE  simulations  over  a 
300  Hz  bandwidth,  centered  at  210  Hz,  and  shaped  by  the 
bulb-signal  spectrum.  As  with  the  real  data,  the  simulated 
data  show  a  direct  arrival  propagating  up  the  array,  a  surface 
path  down  the  array  and  a  bottom  path  traveling  up  the  array 
but  without  the  degree  of  granularity  inherent  in  such  obser¬ 
vations  based  on  large  depth  separation  in  the  VLA.  The  spa- 
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tial  structure  of  the  reactive  intensity  provides  information 
on  the  interference  field. 

The  direction  of  the  reactive  intensity  measured  along 
the  array  gives  an  indication  of  the  divergence  of  Q.  In  a 
region  about  a  maxima  in  the  \p\ 2  field  reactive  intensity  is 
divergent  (V  •  Q  >  0)  and  Q  points  outward  from  the  max¬ 
ima;  maxima  are  sources  of  reactive  intensity.  At  minima  in 
the  |p  | 2  field,  or  nulls,  reactive  intensity  is  convergent 
(V  Q  <  0)  and  Q  points  toward  the  minima;  minima  are 
sinks  of  reactive  intensity.  In  this  coordinate  system,  a 
change  of  sign  in  the  reactive  intensity  component  along  the 
array  from  positive  to  negative  with  increasing  depth  indi¬ 
cates  a  sink  of  Q  while  a  change  from  negative  to  positive 
indicates  a  source  of  Q.  The  data  suggests  a  sink  between 
29-33  m  and  a  source  between  35-39  m,  features  created  by 
the  constructive  and  destructive  interference  of  the  surface 
and  bottom  arrivals.  These  features  are  also  observed  in  the 
simulation.  Sinks  and  sources  of  reactive  intensity  are  defin¬ 
ing  characteristics  of  intensity  structures  in  the  waveguide, 
which  are  further  discussed  in  Sec.  IV. 

Note  that  one  possible  source  of  reactive  intensity  for  a 
time  separated  surface  path  occurs  when  there  is  substantial 
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forward  scattering  from  a  rough  sea  surface.  The  wave  con¬ 
ditions  and  the  propagation  angles  put  the  Rayleigh  parame¬ 
ter,  x  =  2&//rms  sin  9S,  where  0S  is  the  angle  of  the  specular 
reflection  relative  to  the  horizontal,  equal  to  0.36  at  the  high¬ 
est  frequency  (360  Hz).  Although  x  is  less  than  unity,  we 
expect  some  degree  of  forward  scattering  that  could  interfere 
with  the  reflected  surface  path  and  generate  reactive 
intensity. 

B.  Time  averaged  complex  intensity 

At  ranges  greater  than  a  few  waveguide  depths,  a  signal 
requires  a  large  bandwidth  to  exhibit  time  separation 
between  multipaths.  The  frequency  span  of  the  filter  applied 
to  the  experimental  measurements  should  be  set  to  minimize 
measurement  errors;  however,  reducing  the  filter  bandwidth 
has  the  consequence  of  reduced  time  resolution  and  multi¬ 
paths  are  no  longer  resolved.  As  the  bandwidth  is  reduced  to 
form  an  ever  more  narrowband  signal  at  center  frequency  /0, 
the  complex  intensity  loses  its  time  dependence  and  I  and  Q 
can  represented  by  a  single  steady  state  value  (of  positive  or 
negative  sign)  transforming  Ic  into  a  time-independent  func¬ 
tion  of  frequency. 

It  is  worth  further  clarifying  the  role  of  bandwidth  B  in 
constructing  the  complex  intensity  envelope,  and  subsequent 
time  averaging.  Figure  4  illustrates  a  partitioning  of  the 
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FIG.  4.  (a)-(c)  Time  series  of  the  active  (black)  and  reactive  (gray)  intensity 
envelopes  at  32.5  m  depth  for  three  filters  with  B  =  10  Hz  at  three  different 
center  frequencies  for  the  A1  bulb  signal  (range  190  m).  (d)  Time-averages 
of  narrow-band  filtered  signals  at  32.5  m  depth,  B  =  10  Hz,  spanning  a 
300  Hz  band  centered  at  210  Hz.  The  intensity  quantities  here  are  normal¬ 
ized  by  the  maximum  value  of  |p  1 2/po  c  at  the  corresponding  frequency. 


300  Hz  bandwidth  data  previously  shown  in  Figs.  3(a)-3(b) 
into  frequency  bins  located  at/0  by  applying  filters  with  a 
smaller  B  to  the  signal.  Figures  4(a)^l(c)  show  examples  of 
the  data  filtered  with  B  =  10  Hz  at  three  different  center  fre¬ 
quencies,  /0.  Estimates  of  the  component  of  (I)  and  (Q) 
measured  along  the  array  from  a  filter  ( B  =  10  Hz)  are  highly 
dependent  on  the  center  frequency.  For  example,  at 
/o  =  205  Hz,  shown  in  Fig.  4(b),  the  component  of  (I)  ~  0 
and  ( Q )  ~  —0.05  along  the  array  while  at/0  =  255  Hz,  shown 
in  Fig.  4(c),  (I)  «  0.03  and  (Q)  «  0.  Figure  4(d)  shows  a 
summary  of  (I)  and  (Q)  for  the  frequency  bins  spanning  the 
300  Hz  bandwidth.  Although  the  two  processing  schemes 
result  in  very  different  times  series,  the  sum  of  (I)  and  (Q) 
over  the  300  Hz  band  shown  in  Fig.  4(d)  gives  the  same  esti¬ 
mate  of  (I)  and  ( Q )  from  a  single  300  Hz  bandwidth  filter. 

The  time-average  of  the  active  intensity  envelope,  (/), 
represents  the  mean  flow  of  energy  and  is  equivalent  to  the 
time-average  of  instantaneous  intensity,  (S).1  The  time- 
average  of  the  reactive  intensity  envelope,  (Q),  indicates  the 
presence  of  a  non-propagating,  oscillatory  component  of  S 
over  the  averaging  period.  This  oscillatory  component  (e.g., 
here  due  to  multipath  interference)  depends  on  the  degree  to 
which  p  and  v  were  not  in  phase  during  the  averaging  period. 
The  time-independent  Ic ,  as  constructed  from  narrow-band 
(I)  and  (Q) ,  is  used  hereafter  to  describe  the  spatial  structure 
of  intensity. 


IV.  SPATIAL  STRUCTURE  OF  INTENSITY  IN  AN 
OCEAN  WAVEGUIDE 

Acoustic  intensity  in  a  waveguide  will  exhibit  a  spatial 
structure  corresponding  to  regions  of  constructive  and  de¬ 
structive  interference.  These  structures  are  characterized  by 
the  properties  of  the  time-independent  complex  intensity  and 
fall  into  two  categories:  crests  and  vortices. 

Crests  occur  in  regions  of  constructive  interference 
where  there  are  maxima  in  the  \p\ 2  field.  In  these  regions, 
mean  active  intensity  follows  streamlines  propagating 
energy  away  from  the  source;  a  crest  follows  the  most  ener¬ 
getic  streamline  of  J.  Constructive  interference  produces  gra¬ 
dients  in  the  \p\ 2  field,  and  as  V|/?|2  cc  Q,  the  reactive 
intensity  vector  points  outward  from  a  particular  maximum 
in  the  \p\2  field. 

Vortices  occur  in  regions  of  destructive  interference 
where  there  are  minima  in  the  \p\2  field.  The  mean  active  in¬ 
tensity  forms  a  vortex,  the  center  of  which  is  a  null,  or  dislo¬ 
cation,  where  both  \p\2  and  Ic  are  zero.  In  the  vortex  region, 
active  intensity  can  have  a  left  or  right-handed  curl  (V  x  J) 
while  reactive  intensity  converges  toward  the  dislocation. 
Some  experimental  and  simulated  measurements  of  disloca¬ 
tions  are  described  in  D’ spain  et  al4 

Theoretical  work  by  Eliseevnin  et  al.2  shows  the  spatial 
structure  of  active  intensity  for  two  propagating  modes,  relat¬ 
ing  the  spatial  complexity  and  strength  of  the  intensity  vorti¬ 
ces  to  the  difference  in  mode  number.  Although  the  two  mode 
case  is  useful  to  visualize  intensity  structures  like  the  vortex 
region,  the  spatial  structure  of  intensity  due  to  two  modes  is 
not  preserved  when  considering  the  complete  modal  sum. 
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Additional  modes  change  the  interference  pattern,  e.g.,  chang¬ 
ing  the  number  and  location  of  nulls  in  the  intensity  field. 

A.  Spatial  structure  of  intensity  in  a  Pekeris 
waveguide 

Some  additional  insight  is  gained  if  we  examine  spa¬ 
tial  structure  of  intensity  in  the  context  of  a  Pekeris  wave¬ 
guide.  The  complex  intensity  field  of  a  harmonic  source  in 
a  Pekeris  waveguide  can  be  approximated  by  the  sum  of 
trapped  pressure  and  particle  velocity  modes.  The  pressure 
field  is 

N 

p  =  An^n  (,z)Hq^  (knr )  e~iwt ,  (12) 

n=  1 


and  the  horizontal  and  vertical  particle  velocities  are 


where  A Kmn  =  km  —  kn  and  the  pressure  modal  amplitude  Bm  is 


Bn 


2AA(z) 

U\/km 


(19) 


Simplifying  for  just  two  modes,  it  is  seen  that  Eq.  (18)  will 
go  to  zero  (a  pressure  null  or  pressure  dislocation)  at  a  range 
and  depth  combination,  such  that  cos[AKmnr]  =  —  1  when  Bm 
and  Bn  are  of  equal  sign  and  magnitude,  or  cos  [A  Kmnr\  =  1 
when  Bm  and  Bn  are  of  equal  magnitude  and  opposite  sign.  A 
null  at  one  particular  depth  repeats  in  range  at  the  interfer¬ 
ence  wavelength,  2n/AKmn. 

The  active  vertical  intensity  is  given  by 


1  N  N 

I z  =  ^  ^  ^  ^  BmCn  sin[A/cm72r] , 

m= 1 n= 1 


(20) 


Vr 


V— k^An'¥n(z)Hi?\knr)e-iat,  (13) 

f^ilc°Po 


where  An  is  the  modal  amplitude,  kn  is  horizontal  wavenum¬ 
ber  for  the  nth  mode,  and  and  H j1}  are  the  zeroth  and 
first  order  Hankel  functions  of  the  first  kind.  The  mode 
shapes  for  the  pressure  and  horizontal  particle  velocity  (XP„) 
are  given  by 


sin 


k2  -  kl 


(15) 


and  the  mode  shape  for  the  vertical  particle  velocity  (<P  J  by 


where  the  vertical  velocity  modal  amplitude  Cn  is, 


r 


2AwOw  (z) 
7l(JOp0y/k^ 


(21) 


At  the  range  of  the  pressure  null,  7Z,  again  upon  simplifying 
to  just  two  modes,  is  zero  over  the  entire  depth  of  the  water 
column  due  to  the  sine  term. 

Normalized  forms  of  the  complex  acoustic  intensity, 
referred  to  as  the  pressure  intensity  indices,  which  were  first 
introduced  by  Jacobsen  in  the  context  of  airborne  acoustics,2 
are  applied  here  to  underwater  waveguides  to  assess  the  na¬ 
ture  of  the  intensity  vector  field.  The  pressure  intensity  index 
mirrors  the  complex  intensity,  with  the  active  pressure  inten¬ 
sity  index  defined  as 


®„(z) 


—  £2  cos 


Z\f^-k2n 


(16) 


r 


/ 

\p\2/Poc 


(22) 


where  the  origin  of  depth  coordinate  z  is  at  the  surface.  The 
modal  amplitude13  is 


A 


n 


i2nP  o*P„(zo) 


ya(2A) 

2 y/k2 -hi 


+ 


P0 


(17) 


where  h  is  the  depth  of  the  waver  column,  z0  is  the  source 
depth,  pb  and  cb  are  the  density  and  sound  speed  of  the  sedi¬ 
ment,  and  P0  is  the  source  pressure  amplitude  defined  such 
that  An  has  the  dimension  of  pressure. 

The  vector  components  of  Ic  are  found  by  inserting  Eqs. 
(12)— (14)  into  Eq.  (6).  Note  that  the  correct  sign  of  Q 
depends  on  the  sign  convention  of  the  time  dependent  term, 
requiring  the  particular  conjugation  of  Eq.  (6)  for  the  ( e~lcot ) 
convention  used  here.14 

With  the  Hankel  function  written  in  its  asymptotic  form 
(valid  when  knr  1),  \p\2  is  given  by 


1  N  N 
;££*»*  „cos[A  Kmnr\, 


\P\2 


(18) 


n=  1  n=  1 


and  the  reactive  counterpart  defined  as 

Q*  =  —^ — • 

\p\  /Poc 


(23) 


(The  asterisk  is  used  here  to  indicate  that  these  are  normal¬ 
ized  intensity  quantities,  and  does  not  imply  a  conjugation.) 
At  a  pressure  null,  the  characteristics  of  the  vector  quantities 
I*  and  Q *  are  readily  manifested  by  their  vertical  compo¬ 
nents,  /*  and  <2*.  It  is  of  interest  to  examine  the  behavior  of 
/*  as  range  r  approaches  such  a  null  location.  Continuing 
with  just  the  two-mode  case,  expanding  Eq.  (18)  and 
Eq.  (20)  about  a  null  range  point  r0  shows  that  their  ratio 
goes  as  ±  1 1  dr  where  dr  =  r  —  r0  and  the  sign  indicates  the 
sense  of  the  curl  of  the  active  intensity  vortex,  the  center  of 
which  is  the  null. 

In  terms  of  Q*,  we  first  note  that  at  a  pressure  null,  Q* 
is  also  zero,  and  in  this  case  the  more  interesting  property 
is  how  the  index  behaves  as  a  function  of  z  going  towards 
the  null.  It  is  worth  reiterating  that  g*  is  proportional  to 
the  vertical  gradient  of  |p|2.  Therefore,  a  similar  expansion 
of  Q*  about  the  null  depth  point  shows  that  the  ratio  goes 
as  l/dz.  This  indicates  convergence  of  Q  at  a  null. 
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FIG.  5.  (a)  The  narrowband  \p\2  field  in  dB  at  frequency  210  Hz  for  a  unit 
point  source  at  37  m  in  a  Pekeris  waveguide,  (b)  Expanded  view  of  I*  and 
(c)  Q*  for  the  white  outlined  box  in  (a),  (e)  /*  and  (f)  Q*  plotted  as  a  func¬ 
tion  of  depth  and  frequency  at  a  fixed  range  of  650  m. 


The  behaviors  of  /*  and  <2*  about  a  pressure  null  point 
for  the  case  of  multiple  mode  propagation  have  a  structure 
similar  to  the  two  mode  case.  Figure  5(a)  shows  \p\ 2  contours 
for  a  point  source  (210  Hz,  depth  of  35  m)  in  a  Pekeris 
waveguide,  parameterized  by  depth  62  m,  cw—  1500  m/s, 
pw=  1024  kg/m3,  overlaying  an  semi-infinite  half-space, 
with  cb  =  1600  m/s  and  pb  =  1600  kg/m3.  The  white  box  out¬ 
lines  a  region  containing  dislocations  at  a  range  where  the 
pressure  field  is  sufficiently  approximated  by  six  trapped 
modes. 

Contours  of  /*  within  the  white  box  in  Fig.  5(b)  show 
vortex  regions  where  |7*  |  >  1 .  The  vortex  center  is  indicated 
by  /*  changing  sign,  describing  energy  flux  propagating  in 
opposite  directions  on  either  side  of  the  pressure  null.  Simi¬ 
lar  contours  of  <2*  [Fig.  5(c)]  indicate  convergence  of  reac¬ 
tive  intensity  at  these  same  null  points.  The  box  in  Fig.  5(c) 
shows  such  a  convergence,  indicated  by  <2*  changing  sign 
such  that  Q*  is  directed  inwards  toward  the  center  of  the 
null.  The  behavior  of  /*  along  the  horizontal  axis  and  Q* 
along  the  vertical  axis  at  the  nulls  show  the  same  abrupt 
change  in  sign  characterized  by  the  two  mode  behavior. 
A  key  advantage  of  the  index  approach  is  that  the  pressure 
intensity  indices  are  maximum  at  dislocations  in  the  field, 
and  detection  of  these  regions  (where  |J*|  or  \Q*\  >  1)  is  fea¬ 
sible  with  one  directional  component,  e.g.,  <2*  and  /*  as  can 
be  provided  by  a  vertical  line  array. 


Finally,  Figs.  5(d)-5(e)  show  /*  and  <2*,  respectively,  as 
a  function  of  depth  and  frequency  at  a  fixed  range  (650  m). 
These  form  mirror  images  of  the  range  depth  plots  in  Figs. 
5(b),  5(c). 

This  range  and  frequency  relationship  is  described  by 
the  waveguide  invariant,  /?.  Most  often  applied  to  the  range- 
frequency  relation  in  \p\ 2,  uses  of  the  /?  parameter  include 
detection  of  internal  waves,15  active  sonar,16  passive  range 
estimation,17  and  source  localization.18  The  waveguide 
invariant  has  also  been  used  to  describe  the  range  frequency 
relationship  of  the  active  intensity5  and  reactive  intensity.4 
We  demonstrate  next  that  /?  also  applies  to  the  indices  /*  and 
{?*,  with  the  demonstration  limited  here  to  the  vertical  com¬ 
ponents  /*  and  <2* . 

Contours  of  constant  /*  in  range  and  frequency  are 
described  by 

dr  dP 

H — y-dr  =  0.  (24) 

dco  dr 

Increments  of  frequency  dco  and  range  dr  where  /*  has  the 
same  contour  level  are  found  from 


dm=_  dp/dr_ 
dr  drjdco ’ 

In  a  Pekeris  waveguide,  /*  can  be  formed  from  Eq.  (18)  and 
Eq.  (20)  giving 


N  N 

EES  mCn  sin[A7cwwr] 

r*  m=  1  n=  1  7V 

h  -  PQc-fr~N  _  p°c  2)  ’ 

E  E5m5«  C0S[A  Km„r] 

m= 1  n—  1 


(26) 


where  A f  and  V  define  the  numerator  and  denominator, 
respectively.  The  modal  amplitudes  Bm  and  Cm  are  weak 
functions  of  frequency  and  range  independent.  This  holds  for 
any  normal  mode  description  of  the  acoustic  field  in  a  lim¬ 
ited  bandwidth  for  which  Bm  and  Cm  can  be  approximated  as 
frequency  independent.16  With  this  approximation,  the  right 
side  of  Eq.  (25)  is 


drjdr  pf-yf 

dit/du 


(27) 


Next,  the  ratio  of  the  partial  derivatives  Vr/V(0  is  further 
simplified  following  arguments  given  by  Brekhovskikh  and 
Lysanov;19  a  demonstration  is  provided  by  Quijano  et  al.16 
for  striation  patterns  in  \p\ 2.  These  arguments  also  apply  to 
AfrIAf(l)  which  leads  to  the  same  simplification,  and  Eq.  (25) 
thus  reduces  to  the  familiar  operational  relation  for 
namely, 


P  = 


r  dco 
co  dr 


(28) 


The  same  functional  relationship  also  applies  to  Q* . 
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B.  Intensity  structure  in  the  TAVEX  waveguide 

In  this  section  the  range-frequency  relationship  in  /*  and 
g*  is  examined  in  the  context  of  the  TAVEX  waveguide 
using  a  combination  of  field  data  and  parabolic  equation 
(PE)  simulation.  PE  simulations  of  the  complex  intensity 
field  for  a  point  source  in  the  TAVEX  waveguide  allow  for 
finer  spatial  and  frequency  resolution  of  /*  and  g*. 

For  the  PE  simulations,  we  use  the  average  SSP 
[Fig.  2(a)]  and  an  experimentally  derived  bottom  model.  For 
the  latter,  time-separated  bottom  reflected  paths  corrected 
for  spreading  loss  are  used  to  construct  an  empirical  relation 
for  the  magnitude  of  the  bottom  pressure  reflection  coeffi¬ 
cient  (Fig.  6),  where  spreading  loss  and  arrival  angle  are 
derived  using  ray  theory.  The  representative  bottom  model 
(black,  dashed  line)  is  defined  by  cb/cw=  1.06,  density 
1600  kg/m3  and  ocb  =  0.019  dB/m,  a  (the  absorption  coeffi¬ 
cient)  applies  to  the  simulation  frequency  of  210  Hz  and 
originates  from  a  model  for  averaged  properties  for  sedi¬ 
ments  in  this  region.20  We  note  that  historical  sedimentologi- 
cal  data21  and  a  chirp  sonar  bottom  profile  of  the 
experimental  site  are  consistent  with  the  estimated  bottom 
model.  The  PE  simulations  of  intensity  vector  properties  are 
carried  out  by  application  of  finite  differencing  of  a  fine  spa¬ 
tial  sampling.  The  depth  and  range  computational  steps  are 
set  to  be  less  than  10%  of  a  wavelength  of  the  maximum 
simulated  frequency. 

The  time-independent  Ic  (discussed  in  Sec.  Ill)  is  used 
to  construct  estimates  of  /*  and  <2*  with  the  TAVEX  data 
as  a  function  of  frequency  based  on  B  =  1  Hz  with  /0  rang¬ 
ing  from  110-3 10  Hz.  For  ease  in  notation,  /*  and  g*  are 
used  to  denote  the  measurement  along  the  VLA,  which  is 
influenced  by  the  array  tilt  discussed  in  Sec.  III.  The  meas¬ 
urements  in  fact  represent  a  vector  sum  of  the  two  direc¬ 
tional  components  according  to  the  tilt  of  the  hydrophone 
pair  used  in  the  intensity  formulation.  Although  the  tilt  can 
have  a  substantial  effect  on  approximating  the  true  value  of 
/*,  the  horizontal  component  /*  is  typically  directed  away 
from  the  source.  The  direction  of  the  tilt  amounts  to  a  posi¬ 
tive  or  negative  shift  of  the  approximate  /* .  The  normalized 


FIG.  6.  Estimates  of  the  bottom  reflection  coefficient  (x)  for  two  sound 
speed  ratios  cb/cw=  1.04  (dark  gray)  and  cb/cw=  1.065  (light  gray)  over  a 
range  of  bottom  densities  1400  to  1800  kg/m3.  The  representative  bottom 
model  (black,  dashed)  has  cb/cw=  1.065  and  pb=  1600  kg/m3  and 
ab  =  0.019  dB/m. 
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quantity  g*,  on  the  other  hand,  is  quite  robust  to  errors 
stemming  from  array  tilt  as  the  horizontal  component,  g*, 
is  nearly  zero  for  ranges  greater  than  a  few  wavelengths. 
Error  in  measurement  aside,  the  frequency  range  relation¬ 
ship  described  by  the  waveguide  invariant  holds  for  any 
directional  component  of  /*  and  g* . 

Figures  7(a)  and  7(b)  show  /*  and  g*,  respectively,  as  a 
function  of  frequency  for  the  13  vertical  intensity  receiver 
depths  at  the  source  range  (190  m)  analyzed  in  Sec.  III.  Val¬ 
ues  where  \I*  \  >  1  suggest  a  dislocation  is  nearby,  but  the 
depth  resolution  is  not  adequate  to  resolve  it  precisely.  The 
PE  simulations  of  /*  [Fig.  7(c)]  and  g*  [Fig.  7(d)]  pinpoint  a 
dislocation  at  27  m  depth  and  frequency  187  Hz,  indicated 
by  a  convergence  of  g*  in  depth  when  |g*|  >  1,  in  conjunc¬ 
tion  with  |7*  |  >  1 .  Applying  this  observation  to  the  data  sug¬ 
gests  that  there  is  indeed  a  dislocation  between  26  and  33  m 
depth  at  frequency  187  Hz  as  suggested  by  the  combination 
of  |/;|>1  and  g*  being  convergent,  highlighted  in  the  box 
in  Fig.  7(a)-7(b). 

It  is  next  of  interest  to  compare  the  PE  simulation  as  a 
function  of  depth  and  frequency  at  fixed  range  190  m 
[Figs.  7(c),  7(d)]  with  one  as  a  function  of  depth  and  range  at 
a  center  frequency  210  Hz  [Figs.  7(e),  7(f)].  The  center  fre¬ 
quency  210  Hz  is  used  here  because  it  corresponds  to  a  peak 
in  the  bulb  signal  spectrum.  The  figures  show  that  contours 
of  /*  and  g*  are  clearly  related  in  both  range  and  frequency. 
For  this  PE  simulation,  the  frequency  range  relationship  fol¬ 
lows  that  of  the  Eq.  (28)  with  a  value  of  /}  «  1. 

Figure  8  shows  the  result  of  PE  simulations  of  the  spa¬ 
tial  structure  of  /*  and  g*  at  210  Hz  over  the  entire  range 
relevant  to  the  TAVEX  experiment.  A  single  frequency  sim¬ 
ulation  can  be  reconstructed  by  a  small  set  of  fixed-range 
simulations,  but  at  each  fixed-range,  the  simulation  is  com¬ 
puted  over  a  band  of  frequencies.  The  frequency  band  is 
then  interpreted  as  a  range  band  according  to  Eq.  (28)  using 
P  =  1.  The  single  frequency  /*  in  Fig.  8(a)  is  reconstructed  in 
8(b)  by  multiple  frequency  simulations  at  24  discrete  ranges 
(separated  by  the  black  vertical  lines).  The  single  frequency 
g*  simulation,  plotted  in  Fig.  8(c),  is  reconstructed  in  Fig. 
8(d).  In  other  words,  /*  and  g*  in  Figs.  8(b)-8(d)  are  com¬ 
puted  as  a  function  of  frequency  between  175  and  245  Hz  at 
discrete  ranges  spaced  every  50  m.  The  mapping  from  fre¬ 
quency  to  range  around  a  center  frequency  of  210  Hz,  and 
subsequent  reconstruction  of  /*  and  g*  versus  range,  display 
the  same  intensity  spatial  structure  observed  in  the  original 
narrow  band  simulation  shown  in  Figs.  8(a)  and  8(c). 

The  dominant  features  indicated  by  /*  are  dislocations, 
or  nulls,  indicated  by  abrupt  changes  in  sign  along  “bands” 
of  elevated  /*  (where  |/*|  >1).  These  dislocations  often  rep¬ 
resent  a  vortex  in  active  intensity,  occurring  at  the  point 
when  a  band  narrows  before  undergoing  an  abrupt  change  in 
sign.  Bands  of  /*  trace  out  a  pattern  related  to  the  modal 
structure  of  the  field.  In  a  two-mode  propagation,  disloca¬ 
tions  repeat  at  the  interference  wavelength.  In  the  case  of 
multi-mode  propagation,  the  actual  dislocations  may  not 
repeat  at  the  interference  wavelength,  but  the  bands  of  /* 
have  a  characteristic  structure  associated  with  the  multiple 
interference  wavelengths  as  shown  in  Figs.  8(a),  8(b).  Along 
the  bands  of  /*,  there  is  also  a  seam  of  convergence  in 
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FIG.  7.  (a)  Experimental  estimates  of  I*  as  a  function  of  frequency  and  depth  corresponding  to  the  closest  bulb  range  Al;  (b)  corresponding  estimates  of  Q*. 
Whited  out  regions  correspond  to  a  6  m  gap  in  data  due  to  an  inoperable  hydrophone,  (c)  PE  simulated  estimates  of  I *  as  a  function  of  frequency  and  depth  at 
the  Al  bulb  range  (190m)  and  based  on  a  source  depth  of  37  m;  (d)  corresponding  estimates  of  Q*.  The  white  box  highlights  a  dislocation,  (e)  PE  simulated 
estimates  of  I*  as  function  of  range  and  depth  at  frequency  210  Hz;  (f)  corresponding  estimates  of  Q*. 


reactive  intensity,  indicated  by  <2*  >  1  above  and  <2*  <  1 
below  the  seam,  or  white  above  black  in  the  in  Figs.  8(c), 
8(d). 

There  are  some  features  in  the  narrowband  case 
which  are  not  perfectly  recreated  in  the  mapping.  At  the 
interface  between  two  ranges,  the  curves  are  not  perfectly 
aligned,  which  could  be  attributed  to  a  depth  dependence 
of  the  [I  parameter  caused  by  the  downward  refracting 
SSP.22  In  cases  where  the  sign  of  /*  is  opposite,  the  fre¬ 
quency  dependence  of  the  modal  amplitudes  may  come 
into  play,  which  is  neglected  when  deriving  the  wave¬ 
guide  invariant.4 

A  similar  range  frequency  relationship  of  /*  and  Q*  is 
also  observable  in  the  experimental  data  between  two  nearby 
bulb  ranges.  Figures  9(a)  and  9(b)  show  /*  for  two  sources 
deployed  at  A6  (675  m)  and  A10  (698  m),  respectively.  Fig¬ 
ures  9(c)  and  9(d)  show  the  corresponding  Q*.  The  general 
pattern  of  pressure  intensity  indices,  as  indicated  by  regions 
of  elevated  /*  [Fig.  9(a)]  and  convergent  Q*  [Fig.  9(c)] 
(highlighted  by  the  dotted  lines)  at  range  675  m  are  up- 
shifted  in  frequency  by  ^7  Hz  at  range  698  m.  The  fre¬ 
quency  shift  is  itself  somewhat  frequency  dependent  over 
the  175-245  Hz  range,  but  the  average  of  7  Hz  puts  1. 
Depth  dependence  of  /3  is  difficult  to  quantify  with  this  ex¬ 
perimental  data  set  due  to  uncertainties  in  the  source  location 
and  possible  changes  in  the  relative  receiver  positions  or 
SSP  as  the  two  measurements  occur  and  hour  apart. 
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Continuing  with  /?=  1,  let  us  use  the  TAVEX  data  to 
reconstruct  a  nominal  image  of  /*  and  <2*  as  a  function  of 
range  and  depth.  Using  eight  discrete  ranges  corresponding 
to  bulbs  deployments  A1-A8  and  Eq.  (28),  Figs.  9(e)  and 
9(f)  show  the  reconstruction  of  /*  and  <2*  at  210  Hz.  Note 
that  closer  ranges  require  a  larger  frequency  span,  following 
Eq.  (28),  to  fill  in  the  gap  between  sources;  a  variable  fre¬ 
quency  range  (< dco  <2nx  100  Hz)  is  employed  to  accommo¬ 
date  the  different  range  gaps  dr  associated  with  A1-A8. 
Vertical  lines  in  the  reconstruction  in  Figs.  9(e),  9(f)  identify 
a  portion  of  data  used.  For  example,  the  farthest  range 
between  825  and  950  m  is  based  on  data  from  A8  which  rep¬ 
resents  measurement  range  900  m. 

There  is  a  degree  of  ambiguity  in  the  experimental 
source  depth,  source  range  and  array  shape,  and  the  bottom 
model  used  in  the  PE  simulation  can  be  viewed  as  only  ap¬ 
proximate.  These  points  must  necessarily  add  to  the  fact  that 
the  ranges  corresponding  to  bulb  signals  A1-A8  were  meas¬ 
ured  at  different  times  and  thus  the  water  column  SSP  cer¬ 
tainly  differs  slightly  for  each  range.  With  these  limitations 
in  mind,  the  TAVEX  data  reconstruction  does  share  some 
features  with  the  PE  simulation,  particularly  at  short  ranges. 
In  general,  both  the  data  and  PE  simulation  show  more  bands 
and  seams  at  ranges  less  than  500  m  than  at  ranges  greater 
than  500  m.  This  is  because  the  500  m  range  represents  a 
transition  of  the  field  into  the  trapped  mode  regime.  For 
example,  a  normal  mode  solution  at  frequency  210  Hz  (based 
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FIG.  8.  (a)  Single  frequency  PE  simulation  of  /*  at  210  Hz  and  (b)  the  reconstruction  of  /*  at  210  Hz  from  multiple  frequency  simulations  at  24  discrete  ranges 
(separated  by  dotted  vertical  lines)  using  the  TAVEX  waveguide  model  (see  text).  Dislocations  are  located  at  the  sharp  transitions  from  a  large  positive  value 
(/*  >  1)  to  a  large  negative  value  (/*  <  —1).  (c)  Single  frequency  PE  simulation  of  Q*  at  210  Hz  and  (d)  the  reconstruction  of  Q*  at  210  Hz  from  multiple  fre¬ 
quency  simulations.  The  single  frequency  I *  and  Q*  reconstructions  are  performed  by  a  mapping  of  frequency  to  range  using  the  nominal  shallow  water  wave¬ 
guide  invariant  value  (/?  =  1). 


on  the  same  waveguide  parameters  and  computed  using 
the  program  kraken)  is  described  by  six  trapped  modes  and 
converges  to  PE  simulated  field  for  ranges  greater  than 
~500m. 

Finally,  regarding  the  bottom  model,  we  note  that  small 
changes  in  sediment  sound  speed  (^10  m/s)  change  the  rela¬ 
tive  positions  of  dislocations  significantly,  an  observation 
that  could  be  exploited  for  bottom  inversion.  On  the  other 
hand,  small  changes  in  sediment  density  10  kg/m3)  have 
little  effect  on  /*  and  Q* .  Additional  layering  in  the  bottom 
model  can  change  the  position  and  increase  the  number  of 
dislocations  in  the  waveguide. 

V.  SUMMARY 

Vector  intensity  properties  in  a  shallow  water  wave¬ 
guide  have  been  studied  with  experimental  and  synthetic 
data.  The  experimental  field  data  are  from  acoustic  propaga¬ 
tion  measurements  made  during  an  experiment  that  took 
place  near  the  Korean  ocean  research  station  Ieodo  in  the 
East  China  Sea,  in  waters  62  m  deep.  The  synthetic  data  are 
based  on  parabolic  equation  simulations  using  a  geoacoustic 
model  derived  from  the  experimental  measurements. 


Vector  intensity  is  discussed  in  the  context  of  real¬ 
valued  instantaneous  intensity  S(t)  and  complex- valued  in¬ 
tensity  Ic(t).  Since  the  field  measurements  are  made  with  a 
vertical  line  array  it  is  only  the  vertical  component  of  S(t), 
and  the  real  and  imaginary  parts  of  Ic(t)  corresponding  to  the 
envelope  of  active  and  reactive  intensities,  respectively,  that 
are  described  and  interpreted.  By  envelope,  we  mean  a  low- 
passed  version  of  active  and  reactive  intensity  without  mod¬ 
ulation  associated  with  a  center  frequency.  However,  it  is 
shown  that  the  vertical  components  contain  most  of  the 
waveguide  information.  In  particular,  properties  of  the  wave¬ 
guide  are  revealed  by  examining  the  vertical  component  of 
the  time-dependent  active  and  reactive  intensity  envelopes. 
The  complex  intensity  for  time  separated  multipaths  reveals 
the  sign  of  the  propagating  vertical  energy  flux,  the  degree 
of  multipath  interference,  and  hints  at  SSP-dependent  focus¬ 
ing  associated  with  the  direct  arrival.  Time  averages  of  the 
active  and  reactive  intensity  envelopes  allow  for  comparison 
of  time-domain  estimates  with  a  frequency  domain  represen¬ 
tation  of  complex  intensity. 

Normalized  field  indicators,  referred  to  as  pressure  in¬ 
tensity  indices,  are  described  as  functions  of  frequency.  The 
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FIG.  9.  (a)  /*  versus  depth  for  bulb  A6  (range  675  m);  (b)  corresponding  /*  for  bulb  A10  (range  698  m).  (c)  Q*  versus  depth  for  bulb  A6;  (d)  corresponding  Q* 
for  bulb  A 10.  Regions  where  I*  is  elevated  above  the  neighboring  frequencies  and  where  Q*  is  convergent  are  highlighted  by  dotted  lines.  Whited  out  regions 
correspond  to  a  6m  gap  in  data  due  to  an  inoperable  hydrophone,  (e)  Eight  source-ranges  (corresponding  to  bulb  deployment  A1-A8)  and  the  waveguide 
invariant  (/ 3  =  1)  are  used  to  reconstruct  I*  at  210  Hz  as  a  function  of  range  at  the  13  depths;  (f)  corresponding  Q*z. 


pressure  intensity  indices  /*  and  Q *  indicate  two  types  of  in¬ 
tensity  structures:  crests  (regions  of  constructive  interfer¬ 
ence)  and  dislocations  (nulls  from  destructive  interference). 
Crests  have  associated  with  them  a  divergence  in  reactive  in¬ 
tensity,  whereas  nulls  have  associated  with  them  a  conver¬ 
gence  in  reactive  intensity  and  vortex  structure  in  the  active 
intensity.  Crests  are  indicated  by  |J*|  <  1  and  \Q*\  <  1  and 
nulls  by  |J*  |  >  1  and  \Q*  \  >  1.  The  behavior  of  the  index  vec¬ 
tor  magnitude  at  null  structures  is  describable  by  the  single 
directional  component  |7*  |  >  1  and  |Q!  I  >  1.  The  relative 
location  of  the  null  can  be  described  by  the  sign  of  the  indi¬ 
ces;  a  null  occurs  where  <2*  is  convergent  and  where  /* 
undergoes  a  change  in  sign. 

It  was  further  shown  that  contours  of  constant  /*  and  <2* 
in  range  and  frequency  are  related  by  the  waveguide 
invariant,  which  was  then  used  to  map  PE  simulated  fre¬ 
quency  data  at  fixed  range  to  another  range  for  some  center 
frequency.  This  analysis  was  repeated  with  the  TAVEX 
data;  in  this  case  a  reconstruction  of  /*  and  <2*  as  a  function 
of  depth  over  the  entire  ~  1  km  range  was  made  using  data 
associated  with  eight  discrete  ranges  between  200  m  and 
lkm.  The  reconstruction  has  certain  limitations  based  on 
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uncertainties  in  the  experimental  environment  and  geometry, 
nevertheless  the  data  reconstruction  and  the  PE  simulation 
displayed  some  common  features. 

The  formulation  of  complex  intensity  and  its  normalized 
quantities  provide  a  new  look  at  properties  of  the  sound  field 
that  are  directly  related  to  waveguide  properties.  We  antici¬ 
pate  that  improved  knowledge  of  the  vector  intensity  field  in 
a  waveguide  will  find  applications  in  passive  and  active  so¬ 
nar,  direction-of-arrival  sensing,  and  acoustic  inversion  stud¬ 
ies  involving  waveguide  properties. 
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